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The aechanical properties of a bone are a good indicator of 
the health and condition of that bone, and possibly of the 
skeletal systea as a whole. Along the better correlated 
aechanical properties to bone condition are stiffness 
properties. However, no clinical Method is currently available 
to aeasure such properties noninvasively. 

The long bones of the foreara and leg are the aost 
accessible for aechanical testing. Hence, aany investigators 
have concentrated their efforts on these bones. Various 
approaches have been taken involving either ultrasonics or 
iapedance testing. 

One such iapedance aethod was developed by Thoapson*. 
However, aore developeaent is needed before this aethod is 
suitable for routine use in a clinical setting. Much of that 
needed developeaent work is presented. 

A aatheaatical aodel of the vibrating foreara and leg 
systeas is developed. Briefly, the aodel consists of a unifora. 


linear, visco-elastic, Euler- Bernoulli beaa to represent the 

VI 33 HI O < 

ulna or tibia of the vibrating foreara or leg systea. The skin 

H «E 

and tissue coapressed between the probe and bone is represented 


by a spring in series with the bean. The remaining skin and 
tissue surrounding the bone is represented by a visco-elastic 
foundation vith mass. 

An extensive paranetric study is carried out to determine 
the effect of each parameter of the mathematical aodel on its 
impedance response. Tvo accomplishments are obtained as a res.lt 
of the study. First, an increased understanding of the effects 
of the parameters is gained. Second, many qualitative 
relationships between the parameters and the characteristics of 
the impedance curve are derived. 

i systems identification algorithm is developed, and 
programmed on a digital computer, to determine the parametric 
values of the mathematical model which best simulate the data 
obtained from an impedance test. The algorithm is based on 
minimizing the error function; a function similar in form to 
that of a least-squares method. 

Due to the complexity of the impedance equations of the 
mathematical model, the error function is very nonlinear vith 
respect to its parameters. Consequently, the system of equations 
obtained from a least-squares approach, is virtually impossible 
to solve. Bence, an iterative procedure is developed which 
involves the calculation of a change in each parameter which 
brings that parameter closer to its correct value. To start the 
iteration procedure, an initial guess for each parametric value 
is obtained using the relationships derived in the parametric 
study. 

Data from several groups of impedance tests and experiments 
have been made available through personal communication vith 
imes Besearch Center. Among them are (1) in vitro monkey 


experiments, (2) nonbiological tests, (3) Thoapson's origional 
in vivo, human tests, ana (4) aore recent in vivo monkey tests. 

The in vitro aonkey experiaents involve the aeasareaent of 
iapedance of a aonkey foreara in several stages as the alna is 
being excised. The aatheaatical aodel is shovn to be a good 
representation of the physical systea by using it in its 
appropriate fora to siaulate the whole set of experiaents vith a 
consistent set of paraaetric values. The nonbiological tests 
involve the aeasureaent of iapedance of two systems: a "rigid" 
aass and an aluainoa bean. These "known" systens give an 
indication of the accuracy of the iapedance nethod. The use of 
the computer program is ueaonstrated by applying it to the in 
v ivo huoan and aonkey data. 

Several r< : commendations are given. Additional in vitro 
experiaents are suggested to further understand the support 
conditions of the foreara and leg systems. Improvements to the 
testing procedure are also suggested. 

The impedance testing procedure, vith the recommendations 
taken into account, proaises to be a very useful clinical tool 
for aeasuring mechanical properties of bones. 


* Thompson, 6. A., 1973, "In Vivo Determination of Bone 
Properties from Mechanical Impedance Measurement," abstract in 
Aerospace Medical Association Annual Science Meeting, Las Vegas, 
pp. 133-134. 
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CHAPTER I 
INTRO DOCTION 


A. FORWARD 


li The N eed for Measurecent of Bone P roperties 

Nuaerous recent studies have centered on the noninvasive 
aeasureaent of aechanical properties of bones in vivo. Hany 
different approaches have been taken such as inpedance aethods 
and ultrasonic aethods. Sose of these approaches will be 
discassed in Sections I.D anu I.E. Host of these studies have 
been concerned vith various kinds of stiffness aeasureaents; 
usually either aodulus of elasticity (E) of the aaterial of a 
bone or the bending stiffness (El) of a whole bone. These 
stiffness aeasurenents have aany clinical applications. Aaong 
thea are the detection and the aeasureaent of the degree of 
deterioration resulting froa osteoporosis and other bone 

W 

diseases and the aeasureaent of the degree of fracture healing. 
However, relationships between stiffness aeasnreaents and bone 
disorders Bust be knows, to aakc the stiffness aeasnreaents 
applicable. These relationships and their clinical applications 
will be discussed in Section X.C. Before this discussion, 
however, a brief review of anatoay is appropriate. 


B. INkTOMT 


1. The Skeleton 


The skeleton is the set of hones which forn the internal 
framework of the body. The functions of the bones are given by 
Rowe (1972) as follows: 

1. The outward fora of the human body depends on the 
shape and size of the bones# which are the lain 
supporting structures for other body tissues# 
particularly the auscles. 

2. Some parts of the skeleton protect the vital 
organs; for example the bones of the cranium protect 
the brain and the thoracic cage protects the heart# 
lungs, liver and spleen. 

3. By ■eans of the leverage obtained through the 
articulation of the bones with one another at their 
joints# the muscles are enabled to carry out 
movements# including locomotion. 

4. The calcium contained in the bones not only 
strengthens them against stresses and strains but also 
serves as a reserve from which it may be withdrawn 
into the blood stream should the need arise. 

5. The red marrow contained in cancellous bone is the 
tissue from which red and some of the white blood 
cells are developed. 

The skeletal system must be maintained so that these 
functions can operate. Hany diseases are associated with the 
deterioration of the bones# inducing adverse effects on their 
functions. 

Bone# like other tissues, consists of living cells and non- 
living intercellular substance. However# the intercellular 
substance (or matrix) in bone tissue, unlike other tissues# is 
calcified. Calcium salts impregnate the cement substance of the 
matrix thus giving bone its rigidity. Hany bone diseases result 
in a loss of these calcium salts and hence a loss of bone 
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rigidity. 

There are basically four types of bores, characterized by 
their size and shape: long, shor«., flat and irregular. Many of 
these bones have been studied froa a variety of different points 
of view, in terrs of aonitoring bone integrity. The long bones 
in the liabs of the body, however, are of greatest interest for 
noninvasive aechanical testing. Their accessibility simplifies 
testing procedures and their beanlike fora facilitates 
aathematical aodeling. 

2. Long Bones 

The following four definitions are conventional auong 
anatoeists. The tern are refers to the portion of the upper liab 
between the shoulder and elbow, while the tern forearm refers to 
the portion between the elbow and wrist. The tera thigh refers 
to the portion of the lower liab between the hip and knee, while 
the tern leg refers to the portion between the knee and ankle. 

The bones of the ara and forearc, shown in Figure 1. la, are 
the huaerus, ulna and radius. Hote the closeness of the ulna to 
the outer surface of the foreara. Little or no tissue lies 
between the skin and the nlna over aost of its length. Thus the 
construction of the foreara aakes the ulna conducive to 
noninvasive aechanical testing. 

The bones of the thigh and leg, shown in Figure 1.1b, are 
the feaur, patella (knee cap) , tibia and fibula. The tibia, like 
the ulna, is close to the outer surface and is also suitable for 
noninvasive aechanical testing. 
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C. CLINICAL AP PLICATIONS OF STIFFNESS MEASURING TECHNIQUES 


1a Bong Disease 

"Osteoporosis is the ten used to describe a group of 
diseases of diverse etiology which are characterized by a 
reduction in the Bass of oone per unit volute to a level below 
that required for adequate aechanical support function." (Krane, 
1977) . Osteoporosis results in a loss of bone strength due to 
the loss of bone caterial. Although osteoporosis is a very 
connon metabolic disorder, often associated with other 
disorders, the etiology in most cases is not known. Two of the 
aost coaaon types of osteoporosis are disuse osteoporosis and 
senile osteoporosis. 

Disuse osteoporosis results froa a lack of stress applied 
to a bone. The type and degree of stress applied to a bone 
significantly affects the reaodeling of bone. Reaodeling of the 
bone is the continuous lifelong process of the formation and 
resorption of bone Material. A lack of stress applied to the 
bone can result in a decrease of bone material (i.e., resorption 
will exceed formation). 

f 

Disuse osteoporosis occurs in paralytics and bedridden 
patients with diseases not related to the skeletal system. Many 
studies have been done on the effects of immobility, some as old 
as thirty years, e.g., Deitrick, Hhedon and Sborr (1948)* 

Bone mineral losses have also been found to occur in 
astronauts after an extended period of time in a weightless 


5 


environment. The changes in calcium ace clearest in the 84-day 

Skylab mission, see Hhedon et ml. (1976). Drinary calcium 

excretion was monitored and measurements of bone mineral content 

(BBC) were taken of several bones. Orinary calcium excretion 

increased steadily during the first fev weeks in flight, and 

leveled off at about double the value observed during the 

preflight control period, vith no suggestion of decline toward 

the end of the flight. A maximal loss of 7.9 per cent in BBC was 

observed in the os calcis while the radius and ulna did not 

change measurably. Among the implications expressed by Bhedon et 

al. (1977) is the following: 

Since mineral is lost differentially in greater total 
amounts from trabecular areas of bone, one must 
consider the possibility that in very long space 
flights local area losses of mineral of a degree 
equivalent to osteoporosis, visible by ordinary X-ray 
would take place and that the strength of critical 
bones would be endangered. 

Bence, during longer space flights such as a flight to Bars (1.5 
to 3 years duration) , significant changes are expected to occur 
in the long bones such as the radius or ulna and particularly in 
the weight bearing tibia. 

Wbedon et al. (1977) also points out that "urinary calcium 
inflight increased steadily to a plateau in virtually the sane 
pattern and degree as previously seen in bedrest studies." 
Bence, one would expect that results from such studies are a 
good indication of the effect of weightlessness. Ongoing 
investigations are being conducted to study this effect over 
long periods of restraint (six months or more). See Ionng and 
Tremor (1 978) . 

Senile osteoporosis is an osteon rosis associated with 
aging. Although the exact mechanisms which act to induce this 
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osteoporosis are not known, it is believed to be at least 
partially caused by horaonal iabalances which occnr with age, 
particularly with post-nenopausal changes in woaen. 

Other diseases such as rickets and osteoaalacia also result 
in a decrease in strength in bone. These two disease., are 
associated with a defective aineralization of bone aaterial. 

2. Bone Strength 

Each of the bone diseases discussed above results in a 
decrease in bone strength, the force required to fracture the 
bone. Therefore, a aeasuring technique would be valuable. 
However, bone strength can not be aeasured directly except by 
■ethods which entail destruction of the specimen. Therefore a 
noninvasive aethod for inferring bone strength is needed. If 
correlations can be found between stiffness and bone strength, 
then the stiffness aeasureaent s, eentioned in Section I. A, will 
be very useful, once correlations are established to the ;oint 
that bone strength can be accurately inferred then the stiffoess 
■easurenents can be used to: (1) diagnose bone diseases, (2) 
deteraine the extent of the deterioration caused by the disease, 
(3) prescribe treatment and (4) caution patients to avoid 
activities which will induce dangerous stress levels in their 
bones. 

3. C orre lat ion Studies 

Although bone diseases usually affect all of the bones in 
the skeletal systea, long bones are aore accessible for testing. 
Thus aost of the studies have been concerned with long bones. 

Bather (1967a) (1967b) was aaong the first to correlate bone 
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strength to other aaterial and geometric properties of the bone. 
Be ran siaple bending tests on fresh, excised, human long bones 
and found strong correlations between bone strength and each 
"aeasnreable" quantities as age, aodulus of elasticity and bone 
geoaetry . 

Further correlation studies have been perforaed to relate 
bone strength to bending stiffness of long bones. Borders, 
Petersen and Orne (1977) tested fifty-six excised, fresh, canine 
long bones (ulnae, radii and tibiae) in three and four point 
bending. Jurist and Foltz (1977) tested forty-five excised, 
eabalaed, huan i ulnae in three-point bending. In each case, the 
force versus deflection was recorded while the bone was loaded 
to fracture. Statistical correlations were found between bone 
strength and various aechanical properties of the bones. 

These two independent investigations were parallel althou gh 
the spcciaens used in each were substantially different. Their 
findings and conclusions support one another. In particular, 
very strong correlations were found between bone strength and 
bending stiffness for the noraal bones tested. BBC was also 
■ensured near the center of each bone tested. Both studies 
indicate a substantial correction between BBC and both bone 
strength and bending stiffness. 

Thus, correlations have been well established between bone 
strength and bending stiffness for healthy bones. Further 
correlation studies involving various kinds of diseased bones 
are needed to establish the effect of these diseases. It is 
reasonable to expect that good correlations can be found for 
diseased bones, since they exist for healthy bones, i reliable 
aetbod for aeasuring bending stiffness would then be very useful 
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as a non-invasive indicator of bone strength. 

*L . F ractur e Healing 

Another potential use of stiffness neasnreaents is the 
deternination of the extent of fracture healing. A few recent 
studies hare already been done in this area, inong the first to 
investigate the feasibility of soch an application were Caapbell 
and Jurist (1971). They aade iapedance aeasureaents on an 
excised, intact husan feaur and further aeasureeents on the saae 
bone in various injurious conditions, concluding that aethods of 
this type are indeed feasible. Further stodies %«re carried out 
by Barkey and Jurist (1974) and Hoekseaa and Jurist (1977) in 
which resonant frequency was correlated to fracture healing. 

Boargis and Burny (1972) perforaed a theoretical study to show 
the effect of a partially healed section on the aechanical 
response of a bone. Abendschein and Byatt (1972) Bade ultrasonic 
aeasureaents to obtain the aodulus of elasticity of bones in 1 

guinea pigs at various stages in the healing process, thereby 
deaonstrating its variation with healing. 

In aeasuring bone properties for the purpose of aonitoring I 

the healing process of a fracture, it would be advantageous to 
know what the bone properties were before the fracture occurred. 

This, of course, is not possible in a clinical setting. Bowever, 

Borders, Petersen and Orne (1977) found, in the case of healthy 
canine bones that paired bones (right and left bones of the *axe 
type froa one aniaal) have virtually identical aechanical 
properties. If this paired bone relationship holds true for the 
huaan skeleton as well, then aeasureaents taken on a partially 

] 

healed bone can be coapared to corresponding aeasoreients on its 
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paired bone to detsraine the extent of healing. 
Cadaver Ev a luation 



Still another potential use of £jd vivo stiffness 
■easorenents is the skeletal states evaluation of cadavers. 
Hunan cadavers are used quite extensively for inpact safety 
studies. 1 noninvasive screening technique would be very useful 
in deternining the suitability of a cadaver to represent a 
specific population in such a test, llthougb this approach to 
cadaver evaluation is presently not in widespread use, the 
concept was introduced and discussed in detail by Orne (1976). 


I>i OTHERS WORK 


!i Ultrasonics 

It was shown in the last section that bone conditon is 
related to the nechanical properties of the bone. Hany 
investigators have attenpted, with varying degrees of success, 
to neasure these properties in vivo . Two najor types of 
approaches have been taken: ultrasonics and inpedance testing. 

Craven, Costinini, Greenfield and Stern (1973) investigated 
the plausibility of aeasuring the speed of sound in ulnae in 
v ivo using a pulse-echo technique. They showed a significant 
difference in their aeasireaents for bones of two extreae groups 
of subjects: young healthy sales and older (post-aenopausal) 
feaales. Farther investigations using this nethod were carried 
out by Greenfield et al. (1975). They deduced the aodulus of 
elasticity froa the speed of sound, aeasnreaents of geoaetry and 
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bone mineral content. 

Abendschein and Hyatt (1970) measured the longitudinal wave 
speed of standardized specimens of huaan femoral and tibial 
diaphyseal cortices, in this prelininary in vitro study, they 
found correlations betveen wave speed and a fev physical 
properties including aodulus of elasticity. Selle and Jurist 
(1966) aade sinilar aeasureaents on whole excised ulnae and on 
olnae in vivo . The in vivo tests wer^ conducted on osteoporotic, 
diabetic and noraal subjects. Saha and Lakes (1977) investigated 
the effect of the soft tissue on the aeasureaent of wave speed 
in long bones. They concluded that the presence of soft tissue 
has a significant effect on these aeasureaents and therefore 
aust by considered. 

In each of the ultrasonic methods discussed above, 
geoaetrical aeasureaents were required to deduce the aodulus of 
elasticity of the bone being tested. These aeasureaents can be 
very difficult to obtain accurately in vivo and aay even be 
impossible in a clinical setting. I technique for aeasuring 
bending stiffness El such as inpedance testing is a aore 
sensitive indicator of bone condition than aeasureaents of 
either the aodulus of elasticity E, or geometric properties such 
as I since both are usually affected by a bone disorder. 
Purtheraore, bending stiffness was shown in the last section to 
be well correlated with bone strength. 

2. Iap edance 

A variety of experimental procedures and apparatus have 
been used to measure the mechanical impedance of excised long 
bones and intact liabs. Host of those who have attempted to 
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model their system at all have used relatively siaple models 

which do not account for all of the significant characteristics 

of the iapedance curves. Entrekin and Abraas (1976) aeasured the 

mechanical iapedance of the hunan foreara but did not atteapt to 

aodel it. Jurist (1970), Jurist and Kianian (1973) and Speigl 

and Jurist (1975) aeasured the aechanical iapedance of a similar 

systea but used it only as a method of aeasuring the resonant 

« 

frequency which they then related to the aechanical properties 
of the ulna. Doherty, Bovill and iilson (1974) made iapedance- 
like aeasureaents on three excised tibia. They concluded that 
••stiffness K, or dynamic aass a, are more sensitive to changes 
in the physical state of the human long bone than is resonant 
frequency P, due to the functional relationship of these 
parameters", i.e., P is proportional to \ K/H. 

Garner and Blackketter (1975) used a finite eleaent aodel 
to siaulate their iapedance data froa a buaan foreara. This 
procedure involves aany X-rays of the foreara and very careful 
aeasureaent to deteraine its geoaetry. 

Thompson (1973) aeasured the driving-point aechanical 
iapedance of a forearm near the middle of the ulna. He modeled 
it with a fair aaouDt of success as a siaple single-degree-of- 
freedoa oscillator over a frequency range froa 65 to 1000 Hz. 
Thoapson's procedure and apparatus will be discussed further in 
the next section. Orne (1974) presented an improved model 
consisting of a viscoelastic beaa in series with a three- 
paraaeter solid to represent the skin. Orne and Bandke (1975) 
and Thoapson, Orne and Young (1976) improved tie aodel further 
by including a few different kinds of viscoelastic foundations 
with aass to represent the tissue surrounding the ulna. This 
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model has potential but core exanination and modification is 
required before it can be used effectively for clinical 
application. 

E t VIEP&TION TE STS IT A MES RESEAR CH CEHTER 
1 ♦ Ap paratus a nd Procedure 

A noninvasive method for measuring the driving-point 
mechanical impedance* of an in vivo human ulna was developed by 
Thompson (1973) . The sane procedure and apparatus has since been 
Modified and used on nonkey ulnae and tibiae, (Peterson, 1977) . 

The forearn (or leg) is suspended across two aluminum 
supports as shown in Figure 1.2. An aluminum block is placed 
over the wrist (ankle) and secured by two screws. A downward 
force is applied through the humerus (feaur) to hold the 
proxiaal end of the ulna (tibia) in place. 

Specially formed plaster pads were Made by Thonpson for 
each subject he tested. The plaster pads were formed to the 
subject’s wrist and elbow to maximize con fort while Maintaining 
rigidity of the supports. Petersen substituted the plaster pads 
with a firm putty tduct seal) to increase comfort of the 

subject, but with questionable results. 

# 

A filcoxon Research Impedance Head (aodel Z-11) mounted on 
the vibrating shaft of a Ling Altec electro-magnetic shaker is 


* Driving- point Mechanical iapedance is precisely defined in 
Section II. B. Briefly, it is the ratio of the amplitude of the 
force to the amplitude of the velocity of the driving-point of a 
system. 
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applied to the ulna (tibia) through a cylindrical probe. The 
shaker is mounted on one end of a lever vith counter weights 
applied to the opposite end. Various sized weights are used to 
apply and control a constant preload force on the ulna (tibia). 
Preloads ranging froa 200 to 600 gram-force (196x10* to 589x10* 
dyne) are used. 

A schematic diagram of the impedance— measuring system is 
shown in Figure 1.3. A sinusoidal electrical input signal is 
generated by an audio oscillator and fed through an audio 
amplifier to the electro-magnetic shaker. The shaker, which 
works on the same principal as a loud speaker, converts the 
electrical signal to a aechanical vibration of the impedance 
head and probe. The probe, when placed against a forearm (leg) , 
forces the ulna (tibia) to vibrate at the frequency at which the 
audio oscillator is set. 

The force and acceleration signals from the impedance head 
are fed through operational aaplifiers and high pass filters to 
a Hewlett-Packard gain-phase meter (model 3565A) . The gain-phase 
aeter displays the gain (in decibels) and the phase (in degrees) 
of the force signal, in digital fora, using the acceleration 
signal as a reference. Traces of the force and acceleration 
signals are also displayed on an oscilloscope. 

The forcing freguency and the two readings from the gain- 
phase aeter are recorded by the operator at many different 
frequencies over a specified frequency range. Thoapson made 
measurements in the range froa 65 to 1000 Bz. Later measurements 
were taken in the range froa 100 to 3000 Bz. 






2 . 


g Raw Data 


The gain reading froa the gain-phase aeter is in units of 
decibels. A gain aeasareaent in bels is defined as the coaaon 
logaritha of the ratio of the power P, of the electrical signal 
being aeasared, to the power P e , of a reference signal. 
Therefore in decibels, the gain is 

G = 10 log P/P 0 (1.1) 

Since for a given resistance, power is proportional to the 
square of the voltage 

G = 10 log ?*/?! = 20 log V/F p (1.2) 

where V is the voltage of the signal being aeasared and V e is 
the voltage of the reference signal. The gain reading froa the 
gain-phase aeter is the gain of the force signal relative to the 
acceleration signal. Since the force and acceleration are each 
proportional to their respective signals, the gain reading is 
G = 20 log cP/c e a = 20 log F/a - 20 log c/c D (1.3) 

where c and c 0 are the constants of proportionality and F and a 
are the force and acceleration anplitudes, respectively. The 
quantity, -20 log c/c e , is not known. Therefore the iapedance- 
■easuring systea aust be calibrated in order to convert the gain 
reading to an iapedance. 

A snail calibration aass is attached to the iapedance head 

in place of the probe. A gain reading for the aass is taken at 

✓ 

100 Hz. This reading should be independent of frequency (at 
least for relatively low frequencies) since F/a in this case is 
the aass a, a constant. Equation (1.3) applied to the 
-calibration aass is 

G/» * 20 log a - 20 log c/c 0 (l.b) 

where G„ is the gain reading for the aass. The result of 


15 


subtracting equation (1.4) from equation (1.3) is 

G - G m ■ 20 log F/a - 20 log a (1.5) 

Solve equation (1.5) for F/a 

F/a * a aatilog (G-G„)/20 (1.6) 

Equation (1.6) is the ratio of the amplitude of force to the 
aaplitude of the acceleration. Impedance, however, is the ratio 
of the aaplitude of the force to the amplitude of the velocity. 
Since the input force (and hence the notion, if the system is 
linear) is harmonic, the relationship between the velocity and 
acceleration amplitudes is 

a * vp (1.7) 

where p is the forcing frequency. Therefore the impedance is 

Z * F/v * ap antilog (G-G m )/20 (1.8) 

A computer program was developed by Thompson to carry out 
the above computations. The calibration Bass and its gain 
reading are entered into the computer followed by each test 
frequency and its corresponding gain-phase readings. The gain 
reading at each frequency is converted to an iapedance using 
equation (1.8). The phase reading at each frequency is adjusted 
by 900 to account for the difference between the acceleration 
and velocity, i.e. , 

phase of impedance = phase of F/a ♦ 90° 

Finally, the resalts are tabulated and plotted, e.g.» see 
Figures 1.4 and 1.5. 


It the purpose and direction of this w ork 


1. Int erpreta t ion of Impedance Measurements 

The mechanical impedance response of a given system 
contains information about the mechanical properties of that 
system. Bence, Thompson's impedance measuring technique 
described in the last section is potentially a very powerful 
clinical tool for determining bone properties. However, it alone 
is not enough. Thompson's procedure produces an impedance plot 
which must be interpreted to extract the mechanical properties 
of the bone being tested. Two major concepts must be developed: 
(1) an appropriate mathematical model and (2) a systems 
identification technique. 

A mathematical model which accounts for the predominant 
characteristics of the system must be developed. Expressions for 
the mechanical impedance of the model must be derived and 
studied in detail to gain an understanding of its behavior. 
Several versions of the model must be considered to determine 
the importance of each of its paraaeters. 

A systems identification technique must be developed to 
determine the values of the paraaeters in the mathematical model 
for any given test, ihen the values are correct, the model will 

t 

generate an impedance plot which matches the impedance plot of 
the system (i.e., the data from the test) over the frequency 
range of the test. The technique must uniquely determine that 
set of values. Furthermore, it must be systematic enough to 
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program on a digital computer, k user oriented program trill be 
written to eliminate the need for a trained operator. 

k set of Titro impedance tests will be discussed and 
analysed using the systems identification technique. These tests 
will establish some verification of the modeling. 

The ultimate goal* of course, is to achieve a working 
scheme to determine bone properties. The scheme will be applied 
to sets of data from several impedance tests to show how it 
works. 
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CHAPTER II 
BATHEKATICAL MODELS 


Lt THE WEED FOR HATHEHATICAL BODELS 

1 . Constr uctio n and Appl ication 

Beal physical systeas can be extremely conplicated and 
difficult to study. It is therefore advantageous to Bake sone 
siaplifying assusptions about the systea to be studied which are 
approxiaa tely correct, thereby constructing a aodel which 
represents the systea. The aodel can then be studied to gain an 
understanding of the systea. Useful relationships between parts 
of the systea can be discovered as an outcoae of the aodel 
studies. 

It is often of interest to aake a specific aeasureaent on a 
part of the systea being studied. Unfortunately, however, many 
physical systeas, especially biological systeas, cannot be 

f 

disassembled to aake that aeasureaent without destroying the 
systea. Therefore, if a reasonable aodel of the systea can be 
constructed with sufficient correlations established between it 
and the systea, then noninvasive aeasureaents can be Bade on the 
systea which infer the aeasureaent of interest throogh the 
aodel. 
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The accarac; of the assumptions Bade it constructing the 
aodel has significant effects on both the outcone of the aodel 
studies and the accuracy to vhich a measurement can be inferred. 
Therefore, these assumptions should be accurate to construct a 
reasonable aodel. 

The aeasurement of interest here is the stiffness of a long 
bone. The noninvasive measurement being modeled is the 
mechanical impedance vhich vill be defined more precisely in the 
next section. A mathematical model of the forearm and leg vill 
be derived, studied and applied to the measurement of bone 
stiffness in the chapters that follow. 

B. IMPEDANCE 


1 . Definition 

In general, impedance is the ratio of input to output of a 
linear systea. A linear system is one in vhich the output is in 
the same proportion to the input, regardless of the amplitude of 
that input. Bence, impedance is independent of amplitude. If the 
input to a linear system is harmonic, then the output vill also 
be harmonic, possibly vith some phase shift. Impedance then, is 

the ratio of the amplitudes of the harmonic input and haraonic 

* 

output and, mathematically, aust be a complex guantitj to 
account for the phase shift. If the output is taken to be the 
physical response of a specific point in the systea then the 
impedance is said to be the impedance of that point. 
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In a mechanical system, the input is usually a force. 2 The 
corresponding output is the Telocity of the point in the systea 
at which the iapedance is being considered. If the point under 
consideration is the point in the systea at which the force is 
being applied then the iapedance is known as the driving-point 
mechanical iapedance (DPHI). 

For a linear systea, the DPHI is independent of the 
amplitude of the input force. 

2 . Justif ication 

The three types of idealized aechanical elements are: Bass, 
daiper and spring. The behavior of any linear aechanical systea 
can be simulated (over a snail enough frequency range) using one 
or some combination of these elements. Therefore, in order to 
clearly define the behavior of a systea, the three basic 
eleaents must be distinguishable on the response curve of that 
systea in what ever fora it is presented. The response curve can 
be presented in a number of ways. It can be presented as the 
ratio of force to acceleration, velocity or displaceaent. 
Furthermore, it can be plotted on either a linear or a log plot. 

The equation of notion for a force f, applied to each of 
the basic elements is given in Tabel 2.1. If the input force is 
harmonic, then the response will be harmonic, and the following 
relationships bold between the aaplitudes of the acceleration a, 
velocity v, and displaceaent f 

a * p 2 * v * p6 (2.1) 


2 In other types of aechanical systems the input night be, for 
example, a torque or a hydraulic pressure. The corresponding 
outputs in these cases are an angular velocity and a fluid flow 
rate, respectively. 
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vbere p is the forcing frequency. The ratios of the aaplitudes 
of the force to the acceleration, Telocity and displacement are 
easily derivable froa equations (2.1) and the equation of lotion 
for each eleaent. These ratios are also listed in Table 2.1. 

Vote that each of the ratios is proportional to an integer 
pover of the forcing frequency p. Therefore, a log-log plot of 
one of the ratios versus the forcing frequency is a straight 
line. The slope of the straight line is equal to the pover of p. 
For exaaple, the ratio of the force to acceleration for a spring 
is 

P/a * kp-* (2.2) 
Taking the log of equation (2.2) yields 

log P/a « -2 log p ♦ log k (2.3) 
Equation (2.3) is a straight line vith a slope of -2 on a plot 
of log P/a versus log p, i.e., the line sakes an angle of 
arctan (-2) * -63. Qo vith the horizontal. In a siailar manner, 
the slope of the straight line produced by plotting each of the 
other ratios is calculated and listed in in Table 2.1. 

To aaziaize the distinguishability betveen the response of 
the aass, daaper and spring, the response curve aust be 
presented in such a vay to aaxiaize the difference in the slopes 
of the response of each of the three basic eleaents. The slopes 
in each case listed in Table 2.1 reveal that this can be 
acconplished by presenting the response curve in the fora of P/v 
(impedance) rather than P/a or P/6. 
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C.. 155 R ELATIONSH IP BETWEEN THE BATHEHATICAL BODEL AND THE 

PHYSICAL STSTEH 


Id Bodeiing a aecbanical systen, the aodel used Bust, in 
soae sense, resenble the actual physical systea. This 
reseablance Bust be evident to give physical Beaniog to the 
paraaeters of the aodel. The physical paraaeters associated vith 
the aaterial characteristics, as veil as those associated vith 
the geoaetrical characteristics, aust be accounted for in as 
auch detail as the investigator is villing to deal vith. It is 
often appropriate to start vith a aodel vhich accounts for the 
■ost obvious physical paraaeters to gain an understanding of the 
systea, and then to progress to other aodels vhich account for 
soae of the finer details of the physical systea. 

Bodeiing of the foreara systea associated vith the 
iapedance-aeasuring procedure developed by Thoapson (1973) 
(discussed in Section I.Z) vas first atteapted by Orne (1974). 
Orne aodeled the ulna as a onifora, linear, visco-elastic, 
siaply-supported, Euler-Bernoulli bean. The skin vhich is 
coapressed betveen the ulna and the probe vas represented by a 
tri-pnraaeter solid in series vith the bean. The haraonically 
varying load applied by the probe is represented by a 
concentrated force applied to the bean through the tri- 
paraaeter solid as sbovn in Figure 2.1. Orne and Handle (1975) 
iaproved this aodel by including a one-degree-of-freedoa aass 
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with elastic and viscous resistance, uniformly distributed along 
tbe bean to represent the tissue surrounding the ulna. This 
refinement produced the capability of tbe model to account for 
the sub-resonances that are evident in the otherwise siooth 
impedance curves, i further refinement was made by Thompson, 
Orne and Toung (1976) in which the one-degree-of- freedom tissue 
model was replaced by a continuous tissue model, idditional 
refinements involving tbe boundary conditions of the bean will 
be presented here. These models will also be applied, with some 
modification, to tbe leg system as well. 

2s. The Bone 

Several assusptions have been made in modeling the bone as 
a unifora, linear, visco-elastic, simple-supported, Euler- 
Bernoulli beam. First of all, a uniform Euler-Bernoulli beam is 
a beam which is based on the following two assumptions: (1) the 
cross section of the beam does not change along its length, and 
(2) the beam is slender enough that shear deformation is small 
compared to bending deformation. The first assumption is 
obviously not true of bones and will be investigated in detail 
in Chapter IT. The second assumption va^> shown to be true by 
Piriali, Wright and Vagel (1976). The bear is also assumed to be 
linear. This assumption was verified by Thompson when he shoved 
that the DPHI is independent of amplitude of the driving force 
provided that amplitude is small. Finally the beam is assumed to 
be visco-elastic. This is a reasonable assumption since the 
structure of bone material, on the microscopic level, is a 
fluid-filled matrix. 


3. The S upports 


Orne (1974) reasons that the sopports at thv ends of the 
ulna are such that tLa resisting Bonent is negligible and the 
transverse rigidity is aach greater than that of the bone, 
therefore the bone is siaply-supported. However, other aspects 
of the conditions at the supports have not been considered. The 
transverse rigidity of the supports when the plaster pads are 
repl . 2 d by putty is guestionable. A possible aisalignaent 
between the downward force applied through the huaerus with the 
support point of the elbow can conceivably cause an effective 
resistance to rotation at the support. 

Several different classical and non-clasical bean boundary 
conditions are proposed as possibilities for representing the 
notion of the bone at the joints. These include various 
combinations of translational and rotational springs at the 
supports. One special case is considered in which the beam is 
extended past the support to a translational spring to represent 
the possible aisalignaent of the huaerus over the support. 

4. The Skin and Tissue 

It is advantageous at this point to propose two 
definitions. The skin and the thin layer of tissue which are 

coapressed between the bone and the probe will be referred to as 

/ 

the skin . All of the ausculatnre, skin and other tissue 
surrounding the bone will be refered to as the tissue . The lack 
of consistency of these definitions in the literature can be a 
-source of aisinterpretanion. Therefore, the proposed definitions 
will be' used here to insure clarity. 

The tissue aodel is presented by Thoapson, Orne and loung 
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(1976) as "an infinite series of one-dimensional visco-elastic 
rods attached to and vibrating with the nlna and rigidly 
attached to and restrained against motion at their opposite ends 
by the radius." This model is conceptually identical to the 
classical problea of a beaa on an elastic foundation. The 
difference is that the classical foundation includes only a 
stiffness element, vhereas the tissue model includes stiffness, 
damping and mass elements. The tissue model vill often be 
referred to as a visco-elastic foundation vith mass, or simply 
as the "foundation." The shear coupling between adjacent fibers 
of the foundation is neglected. The fixed-end boundary condition 
is replaced by a free-end boundary condition when modeling the 
tibia. 

The skin is represented in Orne’s model by a tri-parameter 
solid, as shown in Figure 2.1. This may seem like a reasonable 
representation since one would expect the skin to exhibit 
damping as well as stiffness characteristics. However, a typical 
set of impedance data from a piece of skin shown in Figure 2.2 
indicates springlike behavior over the entire frequency range. 
(Recall from Section II. B that the DPBI of a spring is a 
straight line with a -45 degree slope.) Therefore the skin vill 
be represented here by a simple spring, as shown in Figure 2.3. 
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D. THE BEHAVIOR OF THE MATHEMATICAL MODEL 


1 . In pe da nee Eq uations and Paraoetric Study 

The aathenatical aodel described in the last section is to 
be studied to gain an understanding of the systea. In order to 
conduct this study, equations for the DPBI of the nodel Bust 
derived. These equations will contain, as one of their 
paraaeters, the quantity to be aeasured, i.e. , the bending 
stiffness of the bone. The equations will be nondiaensionalized 
to reduce the nuober of independent paraaeters and then plotted. 
The nondiaensionalized plots will facilitate the study of the 
■atheaatical nodel. 

These plots can be used to study the nodel in a nuaber of 
ways. They will be used to deteraiue the effects that each of 
the aodel paraaeters have on the plots. Further use of the plots 
will be tore productive if the effects of each paraaeter are 
known. 

Quantitatively, they will be used in generating 
approiiaate, seai-eapirical relationships between the paraaeters 
of the aatheaatical aodel and the characteristics of the DPB1 
plot of that nodel. Relationships of this type will be useful in 
obtaining approxiaations for the values of the paraaeters of the 
systea directly froa its DPBI plot. 

Qualitatively, the plots will be used to aid in deteraining 
which paraaeters to include in the aodel of the systea. This is 
accoaplished by coaparing the DPBI plot of the systea to the 
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■odel plots to distinguish between the paraaeters which are 
essential to obtain an appropriately shaped DPSI plot and those 
which are not. 

The DPHI equations and their plots will be the subjects of 
the next two chapters. 
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CHAPTER III 
I MPEDANCE EQUATIONS 


A. THE GENERAL METHOD FOR DERIVING IMPEDANCE EQUATIONS 


1 . Background 

Driving-point mechanical impedance (DPMI) is the mechanical 
impedance of the point in the system at vhicb the driving force 
is being applied. To derive the DPMI of a mathematical model, 
one nust solve the equations of motion, evaluate the steady- 
state solution for the velocity at the driving-point and take 
the ratio of the force to the velocity. The method for deriving 
the DPMI of the lathenatical aodel described in Section II. C is 
presented in this section. 

Orne (1974) and Orne and Handke (1975) have derived the 
DPMI of a siaply-supported bean on a one-degree-of-freedoa 
visco-elastic-f oundation-with-mass. The analysis presented here 
is aore general in that the boundary conditions are not 
restricted to siaply-supported. Six different sets of bonndary 
conditions are considered; the siaply-supported case and five 
nonclassical cases. A diagram of each case is shown in Figure 
3.1. The visco-elastic-foundation-with-mass is continuous and 
two types of boundary conditions on the foundation are allowed. 
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2i Tke Derivation 

A convenient way to define a coordinate systea on the be a a 
is shown in Figure 3.2. y, (x,t) and y*(z,t) are the deflection 
functions defined for 0<x<a and 0<z<b, respectively, shown 
positive in the figure, where the concentrated force is applied 
at x * a (z = b). The equations of Motion are 

El d*y, /dx* * f[l asy,/dx*dt ♦ p d^y./dtz = p, (x,t) , 0<x<a 

(3.1) 

El d*y z /dz* ♦ y}l ♦ p &2y z /5t* = p z (z,t), 0<z<b 

where 

E is the Modulus of elasticity of the beam Material 
I is the area aoaent of inertia of the cross section 
1 \ is the daaping coefficient of the beaa Material 
p is the Mass per unit length of the beaa 

p, ,p a are the force per unit length of the beam due to the 
reaction of the foundation. These equations are based on the 
visco-elastic uni-axial stress-strain law, i.e. , 9- Ee ♦ V[£. To 

deternine the DPHI, the steady state solutions to equations 

(3.1) are required. These solutions are of the forn 
y, ( x , t) « I, (x) exp ipt 

(3.2) 

y i (z#t) « I; (z) exp ipt 

(i.e., every point in the systea is vibrating at the saae 
frequency) where p is the forcing frequency and I, (x) and I* (z) 
are coaplex amplitudes of the bean vibration. Upon substitution 
of equations (3.2) into equations (3.1), the following ordinary 
differential equations are obtained 

E* I d*I,/dx* - pp*I, * Pi (x), 0<x<a 

(3.3) 

E*I d^r a /dz* - pp*I^ * P*(z), 0<z<b 


where 
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E* ■ E(1 ♦ »lip/E) = E(1 ♦ 2i^p/u)) 
p, (x#t) » P, (r) exp ipt 

p 2 (z,t) » P 2 (z) exp ipt (3.<0 

u) - (tt/L) * -iEI/fi 
; « wl?/2E 

Id the cases where the foundation is not included# 
P, (x) = P 2 (z) = 0. Id cases where the foundation is included 

P» (*) " i (x) 

(3.5) 

P 2 (*) * P?P*Xz(z) 

where p* is the complex# frequency-dependent quantity obtained 
by solving the foundation wave equation 

Ef a*u/£>t* ♦ T)f d^u/d^dt - p f d 2 u/dt 2 « 0 (3.6) 

with the appropriate boundary conditions# as indicated in Pigure 
3.3# where 

E f is the modulus of elasticity of the foundation material 
11; is the damping coefficient of the foundation material 
is the density of the foundation material 
u(?#t) is the displacement function of the foundation 
and the shear stresses in the foundation are neglected. Por the 
fixed foundation 

P* = -p, cotY / Y (3.7) 

for the free foundation 

pf = p; tanr/2 / Y/2 (3.8) 

where 

f = pTT /u>f / 41 ♦ 2i?#p/u>; 

p f is the mass per unit length of the foundation 
uOf is the fundamental frequency of the foundation 
is the damping ratio of the foundation. 

The result of substituting equation (3.5) into equation (3.3) is 
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E*I d*I, /dx* - p*p 2 Y, » 0 

(3.9) 

E*I d*Y 2 /dz« - p*P 2 Y 2 = 0 

where p* = p ♦ p*^ . The solutions to these equations are 
I, (x) « i, sinAx ♦ B, cosax ♦ C, sinhax 4 D, coshax 

(3.10) 

Y 2 (z) * 1 2 sinAz 4 b z cos az 4 C z sinhAz 4 D z cosh A z 
where A* « p*p 2 /E*I 

and , B,, C, , D , , k l0 B it C z and D 2 are eight unknown 

constants which depend on the boundary and matching conditions. 
The deflection, slope, bending nonent and shear force functions 
arc found by using eguations (3.2), (3.10) and the following 
e,(x,t) * ay,/dx e £ (z,t) = oy 2 /dz 

fl, (x,t) = E*I 5 2 y, /dx 2 fl 2 (z,t) = E*I £> 2 y 2 /dz 2 (3.11) 

¥, (x,t) = E*I d*y, /dx 3 ? 2 (z,t) = E*I d*y 2 /dz* 

These functions are evaluated at the point of load application 
(x = a and z = b) and substituted into the following Batching 
conditions 

y, (a,t) 4 j z (b,t) = 0 H,(a,t) 4 8 2 (b,t) * 0 (3.12) 

0, (a,t) - 0 2 (b,t) * 0 ?, (a, t) - V 2 (b, t) * F exp ipt 

These are four of the eight eguations required to solve for the 
eight unknown constants in equations (3.10). The reaaining four 
eguations are obtained by evaluating the appropriate functions 
at x = 0 or z * 0 and substituting thes into the boundary 
conditions listed for each case in Table 3.1. 

For the case where the bean is extended a distance e, past 
the left support, a third deflection function with an additional 
four constants is required on the interval -e<x<0. To deternine 
the twelve constants for this case, an additional four equations 
are required. They are obtained fron the following Batching 
conditions at z * 0 


32 


y. (o,t) « o e, (o,t) - e 3 (o f t) « o 

(3.13) 

13 (0,t) « o a, (0,1) - a 3 (0,t) * o 

The deflection anplitode £, al Ihe point where the load is 
applied is determined by evaluating T ( (x) at x * a or I 2 (z) at 
z * b in eguation (3.10). The DPai of the beam is obtained from 
2* * P/ip S (3.14) 

For the case where a transverse translational spring is in 
series with the beam, the DPHI of the system is given by 

2* » (Z$-» ♦ ip/k)“* (3.15) 

Por the case where the spring is not included in the model, 
Z* * Zg. 

The DPHI associated with each set of boundary conditions in 
Table 3.1 is listed in Appendix A. In each case, the diagrams of 
Pigure 3.1, the boundary conditions of Table 3.1 and the 
eguations of Appendix A are each numbered correspondingly. One 
sample DPHI derivation (case 2) is presented in the following 
section to show how the DPHI eguations of Appendix A have been 
derived from the general method presented in this section. 

B. A SPECIPIC EXAMPLE 


1 . Rotation al Spring on One End 

Case 2 was chosen as an example to demonstrate the method 
used in deriving the DPHI. The support at x * 0 is perfectly 
rigid with respect to translation while the resisting moment is 
proportional to the rotation at that support. The support at 
z * 0 is a simple support, i.e., perfectly rigid with respect to 
translation and no resistance to rotation. These conditions are 


33 


listed in aatheaatical fora in Table 3.1. 

The general solutions to the bean equations (3.1) were 
found in the last section to be given by equations (3.2) and 
(3.10) , i.e. , 

Ji (x,t) * [i, sinAx ♦ B, cosAx ♦ C, sinhAx ♦ D, coshAx] exp ipt 

(3. 16) 

y 2 (Z.t) « [A z sinAz ♦ B 2 cosAx ♦ C 2 sinhAz ♦ D 2 coshAz] exp ipt 
The slope, bending aocent and shear force functions are obtained 
froa the deflection functions (3.16) using equations (3.11). 
These functions are substituted into Batching conditions (3. 12) 
to obtain the following four equations 

A, sinAa ♦ B, cosAa ♦ C, sinhAa ♦ D, coshAa (3.17) 

♦ A z sinAb ♦ B 2 cosAb ♦ C 2 sinhAb ♦ 0 2 coshAb « 0 

A, cosAa - B, sinAa ♦ C, coshAa ♦ D, sinhAa (3. 18) 

- A z cosAb * hi sinAb - C z coshAb * D z sinhAb * 0 

-A, sinAa - B, cosAa ♦ C, sinhAa ♦ D, coshAa (3.19) 

- A z sinAb - B z cosAb ♦ C z sinhAb ♦ D z coshAb * 0 

-A, cosAa ♦ B, sinAa ♦ C, coshAa ♦ D, sinhAa (3.20) 

♦ A 2 cosAb - B z sinAb - C 2 coshAb - D z sinhAb 

« P / E*lA* 

These equations each contain all eight of the unknown constants. 
Vith several algebraic steps, four new equations can be 
generated froo these four equations. Bach new equation contains 
only three of the unknown constants. Add and subtract equations 
(3.17) and (3.19). Add and subtract equations (3.18) and (3.20). 
Divide each of the four results by two to obtain respectively 
C, 6inhAa ♦ D, coshAa ♦ C z sinhAb ♦ D z coshAb « 0 (3.21) 

A/ sinAa ♦ B, cosAa ♦ A z sinAb ♦ B 2 cosAb * 0 (3.22) 
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C, coshla ♦ D, sinhla - C £ coshlb - 0 2 sinhlb 

c F / 2E*IA 3 (3.23) 

1, cosla - B, sinla - 1 2 cos lb ♦ B z sin lb 

- -F / 2E*I l 3 (3.24) 

Hnltiply eguation (3.21) by sinhlb, aultiply eguation (3.23) by 
coshlb and add tbe two results 

C, (sinhla sinblb ♦ coshla cosblb) 

♦ D, (coshla sinblb ♦ sinhla cosblb) (3.25) 

♦ C z (sinh*lb - cosh*lb) « p coshlb / 2E*I1 3 
Recall the following hyperbolic identities 
cosh*B - sinh*B = 1 

cosh A cosh B ♦ sinh A sinh B = cosh(A*B) 

cosh A sinh B ♦ sinh A cosh B * sinh(A+B) 

Boting that a ♦ b * L, eguation (3.25) reduces to 

C, coshlL ♦ D, sinhlL - C z « F coshlb / 2E*IA 3 (3.26) 

In a siailar aanner, aultiply eguation (3.21) by sinhla, 

■ultiply eguation (3.23) by coshla and subtract the second 

result fron the first. Then again using the hyperbolic 

identities given above, the result reduces to 

-C, ♦ Cz cosML ♦ D z sinhAL * -F coshla / 2E*IA 3 (3.27) 

Hultiply eguation (3.22) by sinlb, aultiply eguation (3.24) by 

coslb and subtract the second result froa the first 
A, (sinla sinlb - cosla coslb) 

*■ B f (cosla sinlb ♦ sinla coslb) .(3.28) 

♦ A 2 (sin*lb ♦ cos*Ab) * F coslb / 2E*IA 3 
Becall the following tr igonoaetric identities 
cos*B 4 sin*B « 1 

cos A cos B - sin A sin B ■ cos (A*B) 
cos A sin B ♦ sin A cos B = sin(A«B) 
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Again noting that a ♦ b * L, eguation (3.28) reduces to 

-A, cosU ♦ B, sinAL ♦ A* « F cosAb / 2E*IA» (3.29) 
In a siailar wanner, aultiply equation (3.22) by sinAa, aultiply 
equation (3.24) by cosia and add the two results. Then again 
using the trigonoaetric identities given above, the result 
reduces to - 

A, - A i cosAL ♦ B 2 sinAL * -F cosAa / 2E*IA* (3.30) 
Equations (3.26), (3.27), (3.29) and (3.30), vhich contain only 
three of the unknown constants each, apply *to any bean since 
they have been generated without use of the boundary conditions. 

Substitute the deflection, slope and bending aoaent 
functions into the boundary conditions listed in Table 3.1 for 


case 2 to produce the following four equations 

B, ♦ D, * 0 (3.31) 

-B, ♦ B, = k, (A, ♦ C, ) / E*IA (3.32) 

B : ♦ D, « 0 (3.33) 

-B* ♦ Di « 0 (3.34) 

The equations above are easily solved for B, , D, , B z and 0 2 in 
teras of A, and C ( . The results are 

B, = -k, (A, VC,)/ 2E*IA (3.35) 

D, * k, (A, ♦ C , ) / 2E*I1 (3.36) 

Bi « 0 (3.37) 

D z « 0 (3.38) 


Substitute eguations (3.35), (3.36), (3.37) and (3.38) into 

equations (3.26), (3.27), (3.29) and (3.30) and coabine the 
teras which have the sane unknown constant 

C, (coshAL ♦ k, sinhAL / 2E*IA) (3.39) 

♦ A, k, sinhlL / 2E*IA - C t * F cosbAb / 2E*IA» 

C, = C z coshAL ♦ F cosh la / 2E*IA» (3.40) 
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-A, (cosAL ♦ A, sinAL / 2E*IA) (3.41) 

-C, X, sinAL / 2E*IA ♦ k L « P cosAb / 2E*IA» 

A, ■ A z cosAL - F cosAa / 2E*IA* (3.42) 

Substitute eguations (3.40) and (3.42) into eguations (3.39) and 
(3.41) and again conbine teras which hare the sane unknown 
constant and transfer all known teras to the right hand side of 
the eguations 

k z k, /2E*IA cosAL sinhAL 

♦ C t (cosh*AL ♦ k,/2E*IA sinbAL coshAL - 1) 

= P/2 E* I A* [coshAb - coshAa coshAL (3.43) 

- k, /2E*JA sinhAL (coshAa - cosAa)] 

-k z (cos*AL ♦ k, /2E*IA sinAL cosAL - 1) 

- Ci k,/2E*IA coshAL sinAL 

= F/2E*IA J [cosAb - cosAa cosAL (3.44) 

- k,/2E*lA sinAL (cosAa -coshAa) ] 

The last two sets of substitutions have been carried out in 
such a way to reduce the set of eight eguations and eight 
unknowns to a set of two eguations and two unknowns. 

Again, recall hyperbolic and trigonoaetric identities, but 
this tine in a slightly different fora, i.e., 
cosh* (A*B) -1 * sinh* (1+B) 
cos*(A*B) “1 * -sin* (l*B) 

cosh B -cosh(i + B) cosh 1 * -sinh(A*B) sinh h 
cos B.-cos(A*B) cos 1 * sin(A*B) sin A 
Apply these identities to eguations (3.43) and (3.44) with 
a ♦ b * L to obtain 
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1 2 k, /2E*IA cosAL sinbAL 

♦ C L (siuh* A L ♦ k,/2E*IA sinhAL coshAL) 

= P/2E*IA 3 [-sinhAa sinhAL 

- k, /2E*IA sinbAL (cosbAa - cosAa) ] 
A^sinSAL - k , /2E*IA sinAL cosAL) 

- Cl k,/2E*IA coshAL sinAL 
« F/2E*IA 3 [sinAa sinAL 

- k f /2E*IA sinAL (cosAa -coshAa) ] 


Put equations (3.45) and (3.46) ioto aatrix fora 


[A] (C) - (B) 


where 


m 


k, /2E*I A. 
cosAL sinhAL 


sin*AL 
-k,/2E*IA 
sinAL cosAL 


sinh*AL 
♦ k, /2E*I A 
sinhAL coshAL 

-k, /2E*I A 
cosbAL sinAL 


(C) - 



and 


IB} 


-sinhAa sinhAL 


F/2E*IA 3 < 


-k, /2E*I A 
(cosbAa - cosAa) 

sinAa sinAL 
-k, /2E*IA 
(cosAa - cosbAa) 


L 


> 


Hatrix equation (3.47) can now be solved for 1 2 and C 
Craters rule. The deterainant of aatrix [ 1 ] is 
D * -k ( /2’ J '*lA cosAL sinbAL 
k,/2E*IA cosbAL sinAL 

- (sinh*AL ♦ k,/2E*IA sinbAL cosbAL) 

(sin 2 AL - k,/2E*IA sinAL cosAL) 


(3.45) 


(3.46) 


(3.47) 


2 nsing 


(3.48) 


38 


Multiply oat eguation (3.48) and coabine like terss. The 
determinant thee reduces to 

D «= sinhAL sinAL [k,/2E*IA (3-49) 

(sinhAL cosAL - sinAL coshAL) - sinhAL sinAL] 

The solution to natrix eguation (3.47), with the deterninant J) 
of aatrix [A] defined by eguation (3.49), is 

A 2 = F/2E*I A 3 D 

{-k, /2£*IA coshAL sinAL (3.50) 

[-sinhAa sinhAL - k,/2E*IA sinhAL (coshAa - cosAa) ] 

- (sinb 2 AL ♦ k,/2E*IA sinhAL coshAL) 

[sinAa sinAL - k,/2E*IA sinAL (cosAa - coshAa) ]} 

C z * F/2E*IA*D 

(k,/2E*IA cosAL sinhAL (3*51) 

[sinAa sinAL - k,/2E*IA sinAL (cosAa - coshAa)) 

- (sin 2 AL - k,/2E*IA sinAL cosAL) 

[-sinhAa sinhAL - k, /2E*IA sinhAL (coshAa - cosAa)]) 

The constants A 2 , B 2 , C 2 and D 2 are now known froa eguations 
(3. 50), (3.37), (3.51) and (3.38), respectively. The deflection 

aaplitude $, can be calculated froa either -I, (x=a) or (z=b) . 
Therefore if is used then the constants A, , B,, c, and 0, are 
not needed to calculate £. (The calculation of S using I, has 
been aade as a aeans of checking the following calculations but 
it is not presented here.) 

Substitute eguations (3.50), (3.37), (3.51) and (3.38) into 
the second of eguations (3.10) and evaluate the result at i « b 
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S = F/2E*IA 3 D 

{-k, /2E*IA sinAb coshAL sinAL 
[-sinbAa sinhAL - k, /2E*IA sinhAL (coshAa - cos 3 a) ] 
-sinAb (sinh 2 AL ♦ k,/2E*IA sinhAL coshAL) 

[sinAa sinAL - k,/2E*IA sinAL (cosAa - coshAa) ] 

♦ k, /2E*IA sinhAb cosAL sinhAL (3.52) 

[sinAa sinAL - k,/2E*IA sinAL (cosAa - coshAa)] 

-sinhAb (sin 2 AL - k, /2E*IA sinAL cosAL) 

[-sinhAa sinhAL - k, /2E*IA sinhAL (coshAa - cosAa) ]} 
After several steps of algebra, equation (3.52) reduces to 
6 = F/2 E*IA S D sinhAL sinAL 

([sinbAa ♦ k,/2E*IA (coshAa - cosAa)] 

[sinhAb sinAL - k,/2E*IA 

(sinbAb cosAL - sinAb cosh AL) ] (3.53) 

- [sinAa ♦ k, /2E*IA (coshAa - cosAa) ] 

[sinAb sinhAL - k,/2E*IA 

(sinhAb cosAL - sinAb coshAL) ]) 

Define the following three constants 
ft = k,/2E*IA (coshAa - cosAa) 

p = k,/2E*IA (sinAb coshAL - sinhAb cosAL) (3.54) 

Y = k, /2E*IA (sinAL coshAL - sinhAL cosAL) 

Substitute the expression for the deterainant D, fro* eguation 
(3.49) into eguation (3.53) and replace the appropriate terns 
with ft, p- and Y according to equations (3.54) 

' 6 «= F/2 E* I A 3 

[-(sinhAa ♦ a) (sinhAb sinAL ♦ p) (3.55) 

♦ (sinhAa ♦ a) (sinAb sinhAL ♦ p) ] 

/(sinhAL sinAL ♦ Y) 

Finally, substitute eguation (3.55) into eguation (3.14) to 
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obtain the expression for the DPMI 
Z* * 2E*I* 3 /ip 

{[-(sinhAa ♦ a) (sinhlb sinAL ♦ (3) (3*56) 

♦ (sin^a ♦ a.) (sin^Xb sinh^XL ♦ p) ] 

/(sinh^L siniL ♦ r)}~» 

Equations (3*54) and (3.56) are the expressions given in 
Appendix A for case 2. 


HQN-DI RE NSIOHALIZAT I<. W OF IMPEDANCE EQDATIONS 


1 . Non -di Bens ion ali2ation 




The »ost effective vay of studying the role of each 
parameter in a mathematical nodel is to first nondiaensionalize 
the equations associated with that model, and then perforn the 
parametric study. The set of variables and parameters are 
grouped together in a natural vay to fora a set of 
nondimensional variables and parameters, thereby reducing the 
number of parameters to be studied. 

One very natural and convenient vay to uondimensionalize 
the DPMI of a beam is to fora the ratio Za/K, where Z is the 
magnitude of the DPMI, <x> is the fundamental frequency of a 
uniform simply-supported beam of the same length and K is the 

t 

static stiffness 

R « 48BI/L® (3.57) 
of that same simply-supported beam when centrally loaded. The 
nondimensionalized DPMI will be plotted versos the 


nondimensional frequency ratio p/u) where p is the forcing 
frequency. The nondimensional parameters are listed in Table 3.2 


with their definitions. 

The general forts of the DPHI equation is given in Appendix 

A as 


Z* = 2EI** (1 ♦ 2i£p/o;) / ipf(AL) (3.58) 

where E* has been replaced by E(1 ♦ 2i£p/to) according to 
equation (3.4), and f (AL) is a function of XL involving 
trigononetric and hyperbolic functions and nondiaensional spring 
constants. From equations (3.7), (3.8) and (3.10), XL can be 
written as 

*!• = [(p ♦ p f 9<V)P 2 / EI(1 ♦ 2i*p/u>) ]*A L (3.59) 

where 


gm 


-1/y eoty 
< 2/y tan f/2 

0 


for a fixed foundation 
for a free foundation 
for no foundation 


and 


v = "P/*** / fl ♦ 2i& p/ujf 

A few steps of algebra will produce the following equivalent 
expressions in terms of the nondiaensional parameters 
*L = tt -fp7uj (1 ♦ 2i?p/w)-*A (1 ♦ Hg(T))v»* 

(3.60) 

r = TT P/uj B (1 ♦ 2i^Bp/u))-Vh 

Hultiply equation (3.58) by to and divide by equation (3.57). 
After sone simplification, the result reduces to 
Z ion - -rrn/24 -fp7£ (1 ♦ 2i^p/M */* 

(1 ♦ fl g(f))»* f ■-•{*!) (3.61) 

If the spring in series with the beaa is included, then 
multiplying the impedance equation by to /K will simply change the 
additional tern fro* ip/k to i (pA>) /(k/K) . 

If the boundary conditions of the beam are nonclassical, 
then terns involving spring constants will appear in the 
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function f (iL) . The teras that appear are 
2k/Em® and k/2E*I* 

for translational and rotational springs, respectively (see 
Appendix A ) - In terns of the nondinensional paraneters, these 
teras reduce to 

2k/E*lA» * T / (AL) J (1 ♦ 2i^p/uj) 

(3.62) 

k/2E*IA * B / (AL) (1 ♦ 2i£p/u>) 

For case 6 (see Appendix A) the length of the extended part 
of the beam e, also appears in the function f(AL). However, 
everywhere e appears in the function, 1 also appears. Therefore 
the ratio e/L is taken as the nondinensional paraneter £. 

It is also possible to include danping in the nonclassical 
supports. This is done by adding an inaginary, frequency- 
dependent tern to the appropriate spring constant. Thus k would 
be replaced by k ♦ ipc. In terns of nondinensional paraneters, T 
or R would be replaced by 

T (1 ♦ iC T p/b) or B (1 ♦ iC*p/w>) 

respectively, where the new nondinensional paraneter is 
C T = c t wA or C* * c R u»/k 

In the next chapter, the nondinensionalized DPHI equations 
are plotted for several values of the nondinensional paraneters. 
The plots will be studied and nany relationships between the 
paraneters will be deternined. 
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CHAPTEB IV • 
PARABBTBIC STOP! 


It THE basic sibply-sopported beah 


1 . The Bean 

The bone of a vibrating foreara or leg system is 
represented by a visco-elastic beam. Ideally, this beaa is 
assoaed to be simply-supported. This is an incorrect assumption 
for aany driving- point mechanical impedance (DPHI) tests and 
experiments. However, the siaply-sapported beam will be 
investigated here first and the effect of changing the boundary 
conditions will be deferred to the next section. 

Figure 4.1 is the DPHI plot of snch a beaa vith the driving 
force applied at its center. The corves were generated, allowing 
the bean damping to take on five different values. The 
parametric values used to generate this and all other 
nondinensional plots presented in this chapter, are listed in 
Table 4.1. 

Comparing Figure 4.1 to a typical DPBI data plot shown in 
Figure 1.5, it can be seen that the beaa alone does not produce 
all of the characteristics neccessary to model a vibrating 
forearm or leg system. Other elements must be added to the beam 
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to produce these characteristics. However, it is beneficial to 
study and understand the bean itself before adding on these 
other elements. 

it low frequencies, the carves in Figure 4.1 are 
predominately springlike (i.e., the slope of the carve is 
virtually negative one) with a stiffness egual to the static 
stiffness of the beaa. Thus the aagnitade of the DPHI in this 
region can be approximated by 

= K/Pu>w (4. 1) 
where (PuxvfZ^) is any point on the carve in the low frequency 
range and K, in this case, is 

K = 48EI/L3 (4.2) 

The «iniaan points of the carves appear to occur right at 
the fuDdanental frequency of the bean for all values of the beaa 
damping. The aagnitude of the DPHI at that frequency, however, 
does depend on the beaa daaping. To aid in deteraing the nature 
of that dependence, the concept of an equivalent single-degree- 
of-freedoa oscillator is introduced. 

2 . The Equivalent Sin gle- Degree -of- Freed on Osci 11 at or 

A single-degree-of-freedoa oscillator (SDOFO) is a aodel 
which consists of a aass connected to the "ground" by a linear 
spring and a linear viscous damper as shown in Figure 4.2. Its 
DPHI plot, shown in Figure 4.3, was generated, allowing the 
daaping to take on five different values. 

Dote that Figures 4.1 and 4.3 are identical for frequencies 
alaost an order of aagnitude above their fundaaental frequency. 
Define an "equivalent" SDOFO of a beaa as the SDOFO whose static 
stiffness K, fundaaental frequency to, and daaping ratio £, are 
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equal to those of the beam. Then it can be said that a 
cen trallyloaded sinply-supported bean behaves in the same 
Banner (i.e. # has the saae aagnitude and phase angle of its 
DPHI) as its equivalent SOOFO np to frequencies alaost an order 
of aagnitude above their fundaaental frequency. 

The concept of an equivalent SDOFO is the key to deriving 

sone of the relationships between the paraaeters of the bean and 

« 

the characteristics of its DPMI plot. The relative siaplicity of 
the DPBI equation of a SDOFO facilitates the derivations. 1 
relationship derived between the paraaeters of the SDOFO and the 
characteristics of its DPHI plot will be a good approxinaticn 
for any bean that behaves in a sinilar Banner to its equivalent 
SDOFO in the appropriate frequency range. The relationship aust 
be expressed in terns of K, oj and "$ and these paraaeters aust be 
interpreted properly. One such relationship is the dependance of 
the ainiaua point of the DPHI plot on the bean daaping. Its 
derivation follows. 

The DPHI Of a SDOFO is 

Z* = c ♦ i(np - K/p) (4.3) 

In terns of K, u> and % the DPHI is 

Z* = K/u> [2% ♦ i(p/u>- uj/p) ] (4.4) 

The aagnitude of the DPHI is 

Z = J4$* ♦ ( p/u> - Z>7p )* (4.5) 

To find the frequency at which the DPHI is aininun, take the 
derivative with respect to the forcing frequency p, and set it 
equal to zero 

dZ/dp « K/u>* ( p/u> - iv/p) (1 ♦ ^»*/p*) / 


•{T ? 2 ♦ (p/u> - Zy p) 2 * o 

The only real positive solution equation (4.6) is 


(4.6) 
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P«* * w (4.7) 

i.e., the ainimuQ point of the curve does in fact occur at the 
fundanental frequency regardless of the anount of danping 
present. The aagnitude of the DPHI at that frequency, according 
to equation (4.5), is 

Z M iu * 2 £K/cu (4.8) 

Equations (4.7) and (4.8) hold true for a centrally-loaded 
sinply-supported bean with K interpreted according to equation 
(4.2). 

Hore traditional frequency response curves are given in 
terns of a ratio of deflection £, to static deflection F/K, 
rather than force to velocity, i.e., 

SK/F = 1 / ^ ( 1 - p V^T z '♦ (2~pA,) 2 
For exanple. See Thoapson (1972). In this case, the aaxinua 
point occurs at 

p • k- 

Hence, the frequency at which the naxinun occurs is dependent on 
the danping. It was shown above that the nininun point of a DPHI 
curve of a SDOFO occurs right at the fundamental frequency, 
regardless of the danping. This is an additional advantage of 
presenting the response of a systea as an iapedance. 

3. The Location of the Driving Force 

Figure 4.4 is the DPHI plot of a sinply-supported bean with 
the driving force applied at four different locations along the 
length of the bean. 

Each of the curves have the sane shape np to frequencies of 
at least two tines the fundasental frequency. The upward shift 
in the curves is due to the increase in the static stiffness K 
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of the bean, as the driving force is moved away froa the center. 
One sight expect that equations (4.1), (4.7) and (4.8) are still 
valid in this case provided K is interpreted properly, i. e., 

K = 3EIL/a*b* (4.9) 
h few calculations to coapare these equations to the appropriate 
points on the DPfll plot indicate that they are, indeed, good 
approximations. 

In the high frequency range of Figure 4.4, a second 
resonance appears at about four tiaes the fundamental frequency. 
The centrally-loaded beam does not exhibit such a resonance 
since the anti-symmetric modes of vibration are not excited 
under a symmetric loading. 

B. THE EFFECT OF THE BOONDABT CONDITIOHS 


1 . Qualitative Effects 

Ideally, the bone of a vibrating forearm or leg system is 
assumed to vibrate as a simply-supported beam. 1 discussion 
presented by Orne (1974) indicates that this is in fact true of 
the system involved in the test procedure developed by Thoapson 
(1973) (discussed in Section I.E) . Bovever, subsequent 

modifications to this test procedure may have altered the 

/ 

simply-supported condition of the bone. Therefore, it is 
important to investigate the effect of various boundary 
conditions on the DPHI of a beam. 

Figures 4.5 through 4.9 are the DPBI plots of a beam with 
five different, nonclassical boundary conditions: a rotational 
spring on one end, a rotational spring on each end, a 
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translational spring on one end, a translational spring on each 
end and a translational spring on an extended beaa, 
respectively. In each case, tbe nondiaensional spring constant 
vas allowed to take on five different values while holding the 
danping in the bean and supports at a constant valne. 

A sinple support on the end of a bean has infinite 
resistance to translation and no resistance to rotation. The 
DPMI of a sinply-supported bean was presented in the last 
section. Figure 4.1. 

Adding a rotational spring to a support introduces sone 
resistance to the rotation which can occur at that support. The 
effect on the systen is to stiffen it as indicated by the shift 
upward and to the right of the DPMI curves of Figures 4.5 and 
4.6. 

Adding a translational spring to a support relaxes sone of 
the resistance to the translation which can occur at that 
support. The effect on the systen is to reduce its over all 
stiffness as indicated by the shift downward and to the left of 
the DPMI corves of Figures 4.7 and 4.0. 

Extending the bean past its left support and adding a 
translational spring to its end introduces a non-zero bending 
■onent at the left support. This bending aoaent offers sone 
resistance to the rotation which can occur there Just as does a 
rotational spring. Bence, an expression for an equivalent 
rotational spring (EPS) constant was derived by equating the 
bending aowent at the left .support of the extended bean to the 
■osent caused by the sane rotation applied to the EPS. The 


expression is 

B = 3£*T / (12 ♦ 2t*T) 


(4.10) 
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where the nondiuensionalized paraneters are as follows (see 
Table 3-2) 

B the ERS constant 

T the translational spring constant of the 

spring at the end of the extended bean 
€ length of the bean extension 
The result of solving equation (4.10) for T is 
T * 12R / (3£2 - 2s*R) (4.11) 

The set of four values of B used to generate Figure 4.5 were 
used in equation (4.11) to produce an equivalent set of values 
for T. These values were used to generate Figure 4.9. The DPHI 
curves of Figures 4.5 and 4.9 are virtually identical. 
Therefore, any systea which can be aodeled as an extended bean 
with a translational spring on its end can be aodeled equally 
well as a beaa with a rotational spring on one end provided the 
paraoeters of two aodels are related according to equation 
(4. 10). 

2. R e-no ndiaen sion a lization 

It is apparent that the curves of Pigures 4.5 through 4.8 
are sinilar in chape regardless of the boundary conditions of 
the beaa. The location of each curve on its plot, however, is 
affected by the boundary conditions. To investigate this 
further. Figures 4.10 and 4.11 are generated. Figure 4.10 is 
generated by choosing one curve froa each of Figures 4.1 and 4.5 
through 4.8 and re-nondiaensionalizing it with respect to its 
own static stiffness and fandaaental frequency. (Recall that all 
curves thus far have been nondi aensionalized with respect to the 
static stiffness and fundaaental frequency of a centrally- 
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loaded, siaply-supported bean.) Pigure 4.11 is generated by 
changing the damping value used for Figure 4.10 to a lover 
value. 

Expressions for the static stiffness of a bean with various 
boundary conditions have been derived and are listed in Table 
4.2. The fundaaental frequency in each case is obtained by 
solving the appropriate characteristic equation. The natural 
frequencies of a systea occur when the DPMI goes to zero for the 
case of no daoping. For a bean, this occurs when the function 
f (AL) goes to infinity. Hence, the characteristic equation to be 
solved for each set of boundary conditions is obtained by 
sotting the denoainator of f (AL) equal to zero 3 (See Appendix 
1). The lowest value found for AL is then used in the following 
equation to obtain the fundamental frequency 

u?, = (AL) */L 2 HEI/p = (AL/tr) 2 u> (4.12) 
where uj, is the fundamental frequency of the bean in question 
and w> is the fundaaental frequency of a sinply supported beaa. 

The five curves in each of Figures 4.10 and 4.11 are 
virtually identical up to frequencies of at least two tines the 
fundaaental frequency. Hence, two conclusions can be drawn. 

First, recall that equations (4.1), (4.7) and (4.8) hold 
for a sinple-supported beaa. Then these equations also hold for 
(or are at least very good approximations for) beans with other 
boundary conditons provided K is interpreted according to Table 
4.1 and uj and % are interpreted as fundaaental frequencies and 
daaping ratios of the beans. 

Secondly, the shape of the DPHI curve (in the frequency 


3 Characteristic equations obtained in this way are in agreeaent 
with Gonan (1975). 
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range of interest, i.e., op to frequencies of at least too tines 
the fundaaental frequency of the beaa) is deternined by the 
daaping ratio and the location of the DPBI curve on the plot is 
deternined by the static stiffness and fundaaental frequency of 
tfcr. bean. The stiffness of the boundaries of the bean affect 
each of these three quantities in the sane way as does the 
bending stiffness of the beaa. Therefore the bending stiffness 
and the boundary stiffness have the sane effect on the DPHI plot 
of a bean up to frequencies of at least two tines the 
fundanental frequency of the bean. At very high frequencies, the 
curves begin to deviate fron one another. However, the deviation 
is only significant if the danping is relatively low. Therefore, 
the effects of the bending stiffness and the stiffness of the 
supports of a beanlike structure (i.e., an ulna or a tibia) are 
not easily distinguishable on its DPBI data plot. 

C. THE EFFECT OP TAPER 


1 . Qu a l it ative Effects 

It can be seen fron Figure 1.1 that lcng bones are not 
unifom. Sone long bones, such as the ulna, have very severe 
tapers. It is, therefore, worthwhile to investigate the effect 
of taper od the DPBI plot of a beaa. 

A netbod of coaputing the DPBI of a tapered bean is given 
in Appendix B. This nethod was ^ed to generate DPBI plots for 
beaas with two different types of tapers: a linear taper which 
roughly approxinates an ulna and a quadratic taper which roughly 
approximates a tibia (see Figure ft. 12). 
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The nondimensionalized DPHI plot in each case turned out to 
be identical to that of a uniform beam. Apparently, DPHI data 
only provides information about the overall stiffness of a 
beamlike structure and not about its distribution. Therefore, no 
information concerning the nature of the taper of a bone can be 
extracted from its DPHI plot alone. However, using a model in 
which the bone is assumed to be uniform, the average bending 
stiffness is determined. This is the same average bending 
stiffness which was measured and correlated to breaking strength 
in the investigations by Borders, Petersen and Orne (1977) and 
Jurist and Foltz (1977). These correlations provide a means of 
inferring breaking strength from a measurement of bending 
stiffness. Thus, knovlege of the exact geometry is not needed. 

' D. THE EFFECT OF THE FOUNDATION 


1 ♦ Qualitative Effects 

The tissue surrounding the bone of a vibrating forearm or 
leg system is represented by a visco-elastic foundation with 
mass* The boundary of the foundation is either fixed or free as 
discussed in Section II. C. Figures 4.13 and 4.14 are DPHI plots 

of a sinply-supported bean on a fixed foundation while Figures 

/ 

4.15 and 4.16 are DPHI plots of a simply-supported beam on a 
free foundation. Figures 4.13 and 4.15 were generated with the 
damping in the foundation held constant while allowing the mass 
per unit length of the foundation to take on five different 
values. Figures 4.14 and 4.16, on the other hand, were generated 
with the mass per unit length of the foundation held constant 
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vhile allowing the damping in the foundation to take on four 
different values. In each case, the stiffness of the foundation 
is chosen to produce a fundamental subresonant frequency for the 
foundation of one-half the fundaaental frequency of the beaa. 
The arbitrary factor of one-half sufficiently seperates the 
subresonant frequency of the foundation fron the resonant 
frequency of the bean to distinguish their effects. 

The foundation exhibits two major effects on the DPMI 
curves. First, the daxping in the foundation smooths out the 
DPMI in inch the sane way as the danping in the bean. The 
aininun point of the curve noves upward as damping increases 
regardless of the source of the danping (beam or foundation). 
Secondly, the DPMI curve changes drastically in the region 
around the subresonant frequency. This disturbance in the 
otherwise snooth curve is evident in many of the data sets fron 
DPMI tests. It is therefore essential to include a foundation in 
the aatheaatical aodel. 

Zi. Quan tification of the E ffect on the Minimum Point 

Mote fron Figures 4.13 through 4.16 that the aagnitnde of 
the DPMI at the aininun point of the curves is very dependent on 
both the Bass per unit length |i^ , and the damping ratio of 
the foundation. This dependence, expressed in mathematical fora, 
can be used to determine approximate values for these parameters 
for a forearm of leg system directly from its DPMI data plot. 

The fundamental frequency uj ( , of the foundation also 
affects the minimum point. However, it will be useful later to 
have relationships expressing the dependence of p # and on the 
minimum DPMI while holding u> ( constant. The disturbance which 
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appears in aany data sets occurs at approxinately one half the 
fundamental frequency of the beam. Therefore, the relationship 
to be derived will be based on a frequency ratio uj/aJf, of two. 

Due to the complexity of the DPHI equations, the exact 
expression for the dependence of and on the minimum DPHI 
can not be determined. Therefore, approximate relationships are 
derived. The details of the derivation are given in Ippendix C. 
The relationships expressing the dependence of p f and on the 
■inimuB DPHI are 

Z^io/K = 2t ♦ 0.25 ft* p f /p (4.13) 

Z» tw u)/K = 2% ♦ 0.75 p f /p (4.14) 
for the fixed and free foundation, respectively. Since these 
relationships are approxiaate, it j; beneficial to desonstrate 
their accuracy. This is done in Pigu o 4.17. The ainiaua DPHI's 
tabulated in Table 4.2 are shovn as squares on the plot vhile 
equations (4.13) and (4.14) are shovn as solid lines. The 
approximation is quite accurate for the range of values under 
consideration. 

h. THE EFFECT OF THE SPR INS-IN- SERIES 


1 . Qualitative Effects 

t 

The shin of the vibrating forearm or leg system is 
represented by a transverse spring in series with the beam. DPHI 
plots of a simply-supported beam vith the spring in place are 
given in Figures 4.18 and 4.19. Figure 4.18 vas generated vith 
the damping of the beam held constant vhile allowing the spring 
stiffness to take on five different values. Figure 4.19, on the 
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other band, was generated with the spring stiffness held 
constant while allowing the damping of the bean to take on five 
different values. 

At very low frequency, the curves are predominantely 
springlike (i.e., the slope of the carve is virtually negative 
one) . The apparent stiffness is simply the combined static 
stiffnesses of the beam and spring in series. At very high 
frequency, the curves are again springlike. However, the 
apparent stiffness is higher than the apparent stiffness in the 
low frequency range. In the high frequency range, the beam DPHI 
is predoainantely Basslike (see Figure 4. 1) while the spring, of 
course, is still springlike. Thus, the bean DPHI is auch higher 
than that of the spring. Recall that DPBI*s in series add 
according to 

Z * = (1/S* ♦ 1/ZJ)“» (4.15) 
Tne lover of the two DPH 's, the DPHI of the spring in this 
case, dominates the overall DPHI. Therefore, at very high 
frequency the overall DPHI is siaply the DPHI of the spring. In 
other words, the beaa, due to its Bass and daaping, does not 
vibrate at high frequency. 

The apparent stiffness at low and high frequencies have 
often been used to approriiate the bone and skin stiffnesses of 
forearm or leg systems directly from the DPHI plots. A data 
point is chosen from each of the (low and high) frequency ranges 
and used in the following formulas 


^ — Pm an 

n * (1/ZuwPiXK - i/M-» 


(4.16) 

(4.17) 
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Where k = skin stiffness 

K * bone stiffness (eg., 3EIL/a 2 b 2 
for a sinply-supported beai) 

= a data point froa the high frequency range 

(P uNvrZuM,) * a data point froa the low frequency range 

(see Figure 4.20) 

However, large errors are easily introduced with iaproper 
choices of the d? ta points. Recall that the data points Bust be 

j 

taken fron sections of the data plot where the frequency is low 
enough or high enough to indeed produce a slope which is | 

virtually negative 45 degrees. This stipulation does not present 
a problea in the low frequency range. However, the data froa 
aost DPHI tests have not been taken in a frequency range high 

enough to attain the required negative 45 degree slope. However, 

a new relationship has been discovered which allows the skin 
stiffness to be approxioated using the Maxima point (see Figure 
4.20) which occurs just before the high frequency negative slope 
on the data plot. This eliainates the need for the high 
frequency data. 

2i Q uanti fi cation of the Effect on the Maxiaua Poin t 

It can be seen froa Figures 4.18 and 4.19 that the aaxiauB 
point is severely affected by the spring. Although the Baxiaua 
value of the DPHI has a significant dependence on the daaping of 
the beaa, the frequency at which it occurs does not. Therefore, 
an approxiaate relationship between the stiffness of the spring 

| 

and the frequency at which the aaxiaua DPHI occurs can be 
derived which is independent of the daaping in the beaa. 

To find this relationship, two siaplifications are 




1 1 jadflwtfn rt .lirrr fi ; 
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introduced to facilitate the analysis. First, replace the beam 
by its equivalent single-degree-of-freedoa oscillator (SDOFO) . 
Recall fro® Sections IT. i and IV. B that a bean, regardless of 
its boundary conditions, behaves in the sane Banner as its 
equivalent SDOFO up to frequencies of at least tvo tines their 
fundamental frequency. In aany cases the similarity in behavior 
extends to as high as an order of magnitude above the 
fundamental frequency. Recall further that at high frequency, 
the DPHI of a spring-in-series doninates the total DPMI. 
Therefore, a spring in series with a bean behaves in the sane 
■anner as a spring in series vith a SDOFO at any frequency 
provided the spring is soft enough. 

Secondly, since the frequency of interest is assumed to be 
independent of the daaping, set the damping equal to zero. Then 
the frequency vhich Bakes the DPHI aaxiau® vill actually be the 
frequency vhich makes the DPHI approach infinity. Thus, the 
nodel to be analysed is that vhich is shovn in Figure 4.21 vith 
% = 0 . 

The DPHI * s of the SDOFO and the spring are, respectively 
Z* = nip ♦ K/ip <4. 18) 

If = k/ip (4.19) 

vhere p is the forcing frequency. The overall DPHI, according to 
equation (4.15) is 

Z* = (1/(aip ♦ K/ip) ♦ 1/(k/ip) )"* (4.20) 

After replacing ■ by K/u)* and performing several steps of 
algebra, equation (4.20) becomes 

Z* * -iK/uj (k/K)/(pA>) (1 - pVa*) / (1 ♦ k/K - p*A>*) (4.21) 
The DPHI approaches infinity vben the denominator of equation 
(4.21) approaaches xero 


58 


1 ♦ k/K - p*,/u,2 * 0 (4.22) 

Therefore the frequency at which the DPMI is naxioum is given by 
Pd»A> 2 • 1 ♦ k/K («-23) 

Solve equations (4.17) and (4.23) simultaneously for k and K 

* = ZutwP^w PZo/u* (4.24) 

K * Zu**Pu,w PL,A’ 2 / (P5../u; 8 - 1) (4.25) 

Equations (4.24) and (4.25) can now be used to approximate the 
bone and skin stiffnesses without the use of equation (4.16) # 
i.e., without the use of a data point froa the very high 
frequency range. 

Figures 4.18 and 4.19 show that the location of the ainiaua 
point is only slightly affected by the presence of the spring. 
This indicates that the relationships discussed in Section I?. A 
and IV. B (equations 4.7 and 4.8) which relate the ainiaum point 
to the danping ratio and the fundaoental frequency are still 
approxinately valid in the presence of the spring. This is also 
verified by considering the frequency which Bakes equation 
(4.21) go to zero, i.e., set the numerator equal to zero 

1 ~ P *,„/«;* = 0 (4.26) 

or 

P* M /w* * 1 (4. 27) 

Since equations (4.23) and (4.27) were obtained by 
considering the case where £ = 0, they are approximations which 
are independent of the beaa daaping. To investigate the accuracy 
of these approximations, the ainiBua and aaxiDua points of the 
DPBI of the aodel of Figure 4.21 can be found without setting 
the beaa daaping equal to zero. Although this analysis is nearly 
icpossible in closed fora, the first few teras of a Taylor 
series solution can be found. This very lengthy analysis is 


outlined in Appendix D. The first three terns of the solutions 
are 

ftix * s ♦ 1 ♦ 2/S (2*S) / (14 S) ** (4.28) 

- 2/S 3 (2 ♦S)/(US>) (406S* 13S 2 *4S 3 ) ♦ ... 

ft** * 1 - 4/S 5* ♦ 8/S 3 (2+3S) - ... (4.29) 

where S = k/K and 0 * pA>. Both series converge for 0<£<1 and 
S>1 which is the range of values of interest. 

Several typical values of % and S have been tried in 
equations (4.28) and (4.29) and compared to the results from 
equations (4.23) and (4.27), respectively. For exanple, with $ - 
0.2 and S = 5, equation (4.28) yields * 2.453 while equation 
(4.23) yields = 2.449. This and *any other sets of values 
indicate that equations (4.23) and (4.27) are indeed very good 
approximations. 
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CHAPTER ? 

THE S YSTEMS I DENTIFICATION ALGORITHM 


ii THE NEED FOR A SYSTEMATIC METHOD 


V. The Need 

Each of the parameters of the natheaatical aodel 
corresponds to one (or some coabination) of the geoaetrical or 
■aterial properties of the vibrating forearn or leg systen. The 
driving-point mechanical iapedance (DPMI) of the systen is 
■easured in a vibration test (Section I.E). The DPMI of the 
■odel is calculated and depends on the values chosen for its 
parameters (Section III. A). Therefore, the set of paraaetric 
values for the aodel which generates a DPMI carve that closely 
coincides with the DPMI data points of the systen infers the 
geoaetrical and naterial properties of that syaten. A nethod for 
finding this set of paraaetric values is needed. 

2. Req uirements 

To obtain a consistent interpretation of the DPMI data, the 
method used to find the paraaetric values (hence forth referred 
to as "the nethod") nust Joe repeatable and systematic. The 
nethod nust be repeatable in the sense that each tine it is 


applied to a given set of DPHI data it aust produce the saae 
results. The aethod aust be systematic enough to progran on a 
digital computer for on-line analysis. 

Although computers are capable of perforaing tremendous 
amounts of computation, they are incapable of Baking subjective 
decisions. The aethod aast be completely objective in nature and 
expressible in mathematical fora. 

Pinally, the computer program which employs the aethod must 
be set up in a user-oriented fashion. The user in a cliiical 
situation should not need extensive computer experience in order 
to easily obtain results. 

IL THF e rror f uncti on 

1 . Definition 

The first step in developing the aethod is to define an 
error function which quantifies the difference between the 
measured DPHI data and the calculated DPHI of the mathematical 
aodel. The parametric values of the model will then be chosen in 
a systematic way to minimize the error function. This is 
accomplished using a systems identification algorithm (SIDA) 
which is analogous to the classical least-squares approach to 
curve fitting. 

The error at frequency p n , is the difference between 
the measured DPHI 2„ , and the DPHI calculated using the aodel 
Z n (Pi), as shown in Pigure 5.1a. The error function E, is the 
finite sub over all the discrete test frequencies of the squares 
of the percentage errors e*/Z n , divided by the number of data 
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points. The percentage error is used rather than the error 
itself because of the wide range of absolute values which the 
DPMI can take in a single DPHI test. The division by the number 
of data points normalizes the error function so that a 
comparison of its value from two sets of data with different 
numbers of data points is aeaningful. The error, and hence the 
error function, is a function of the parameters of the model 
since it depends on the DPUI of the model. An example of an 
error function as a function of one of the model parameters, 
represented by P t , is shown in Figure 5.1b. 

2 . Analysis 

(Mathematically, the error function is expressed as 

E = 1/N £ [ (2„ - Zn(Pi))/Zr>]* (5.1) 

«•! 

where N is the number of data points. To obtain the parametric 
values using a classical least-squares approach, one would set 
the derivatives of the error function with respect to each of 
the parameters equal to zero. The resulting equations would then 
be solved directly for the parametric values. Due to the 
complexity of the function which represents the DPHI of the 
model, however, this approach is impractical if not impossible. 

Since the DPMI of the model is a continuous function of the 
model parameters, it can be expanded in a Taylor series. 

e= i/n£i/z* [£„ - (Z 0 ♦Tdz./dPi AP L ) y (5.2) 

where H is the nuaber of model paraseters. Higher order terms of 
the series have been neglected and the function which represents 
the DPKI of the model and its derivatives are evaluated at some 
initial set of estimated parametric values. 

Using this form of the error function, changes in the 
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paraaetric values AP ; , rather than the paraoetric values 
thenselves, can be chosen to ainiaize the error function. To 
accoaplish this, set the derivatives of the error function vith 
respect to the changes in the paraaetric values egual to zero 
d E/ dAPj « -2/N £ 1/Z* [Z n - (Z n ♦ £dZ„/dP ; APJ ] dZ n /dPj 

* 0 ; j = 1,2.. .H (5.3) 

After a fev steps of algebra, equation (5.3) becoaes 

dE/dAPj = - 2/H (£ n - Z„)/Z2 dZ„/dPj (5.4) 

- 12 1/Z* dZ n /dP t dZ n /dPi AP L ] = 0; j = 1,2...H 

n-» 

Therefore the equations to be solved are 

[A] (AP) = (B) (5.5) 

ihere the coaponents of the aatrices are 

Ai_j = X VZ? dZ„/dP l dZ^/dPj (5.6) 

and R v = ^ “ Zn)/Z^ dZ n /dPj (5.7) 

n •>< 

The derivatives of the DPHI of the aodel vith respect to each of 
the paraueters is given in Appendix E. 

3. Application 

Since changes in the parasetric values are calculated 
rather than the paraaetric values theaselves, the procedure is 
iterative. The coaponents of the A and B aatrices are calculated 
using the paraaetric values obtained froa previous iteration. 
The changes in the paraaetric values are calculated froa 

(AP) * [A]-* (B) (5.8) 

and added to the old set of paraaetric values to obtain a new 
set. Each succesive set of paraaetric values will reduce the 
value of the error function. The procedure is repeated as aany 
tiaes as necessary to obtain an acceptable set of paraaetric 
values. 
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To begin the iterations, an initial set of paraaetric 
values nust be chosen which will facilitate quick convergence. 

C. CONVERGENCE AND THE INITIAL GPSSS 


1 . Definition 

Convergence is said to have occurred in an iteration scheae 
when further iterations no longer improve the result. In terns 
of the SIDA, convergence has ocurred when the relative change in 
any given paraaetcr becoies sialler than a specified amount, 
e.g., 0.1 percent. The characteristics of the error function 
have a considerable effect on the convergence of the SIDA. 
Therefore, soae control Bust be maintained over the error 
function to insure convergence for the DPHI data from any 
fcreara or leg vibration test. 

2 . Restrict ions on the Bathesatical Hodel 

When two paraaeters of a given Bathesatical nodel have very 
similar effects on its DPHI curve, the effect of changing one 
paraaetric value Bay cancel an opposite effect in the other to 
produce no net effect in either the DPHI curve or the value of 
the error function. In this case, the error function Bay contain 
an infinite number of Bininua points along soae curve in the 
error function space. There is no way to distinguish between 
these ainiaua points. Therefore the DPHI data does not contain 
enough inforaation itself to uniguely define all of the 
paraaeters of the aodel, and the SIDA will diverge. This problea 
has ocurred with the boundary conditions of the bean and with 
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the daaping. To eliainate the problea, soaething sore Bust be 
known about one of the two paraneters. A constant value can then 
be assigned to it, allowing the rest of the paraneters to be 
deterained by the SIDA. 

It was shown in Section 17. B, tbat the static stiffness of 
a bean, and hence its DPHI, is affected in such the sane way by 
the bending stiffness of the beai itself as by the stiffness of 
the boundaries. Therefore, if the aodel inclndes a spring at one 
or both ends of the bean. Then the DPHI data does not contain 
enough inforaation to deteroine all of the paraoetric values. 
Therefore, the characteristics of the supports of the foreara or 
leg aust be known a priori . One way to avoid the necessity of 
determining the support characteristics is to always plac the 
foreara or leg in the fixture in such a way to insure that the 
supports are virtually simply-supported. 

The sharp peaks of the uiniaua and aaziaua points of the 
DPMI curve of an undaaped beaa are rounded-off when damping is 
added. Tht extent of the rounding-off depends on the aaount of 
daaping present but not on its location, i.e., in the beaa or 
foundation, as was shown in Section IV. D. Since both the bone 
and the tissue contribute to the overall daaping of the systea, 
the DPMI data does not contain enough inforaation to deteraine 
all of the paraaetric values. A constant value will be assigned 
to one cf the daaping ratios, tbns allowing the other to be 
deterained by the SIDA. It will be seen in Section TI.A that the 
tissue contributes aucb aore to the overall daaping than does 
the bone. Consequently, the DPBI is relatively insensitive to 
the value chosen for the daaping ratio of the beaa. Therefore, 
it will be held constant at five percent of critical daaping in 
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the fundamental sole, a reasonable value. 

With the boundary conditions being specified and the 
damping in the bean held constant* the aodel has six parameters 
to be deterained by the SIDA. They are the bending stiffness El* 
and the fundamental frequency to, of the bean; the aass per unit 
length the fundamental frequency ia (t and the damping ratio 
, of the foundation and the stiffness k, of the spring. This 
version of the model will be referred to as the six- parameter 
model (6PM) . 

3 . The Initial Guess 

Even if all of the parameters are such that their effects 
on the DPMI curve are independent* it is possible that more than 
one set (but not and infinite number of sets) of parametric 
values exist which will minimize the error function for the DPMI 
data from any given vibration test. One of these referred to as 
the correc t solution * is the set of parametric values 

corresponding to the true geometric and material properties of 
the forearm or leg system being tested. 

Several successive iterations of the SIDA can produce a set 
of parametric values associated with one of the local minimum o’" 
raximum points of the error function. To illustrate this 
concept* an error function is shown in Figure 5.1b. Orly one ci 
these minimum points represents the correct solution* and it 
appears to be the only one in which all of the parametric values 
are positive. The initial values chosen for the parameters to 
start the iterations, referred to as the initial guess * 

determine whether or not the SIDA will converge and to which 
minimum or maximum point. Therefore, the initial guess must be 


close enough to the correct solution to allow the SIDA to 

converge to it. 

The leans for acquiring the initial guess is provided by 
the relationships established in the paraoetric study (Chapter 
IV) . The initial guess is calculated fro* a few Key data points 
using these relationships. In lany cases, the intial guess is 
close enough to the correct solution. However, if one or lore of 
the Key data points happens to contain an excessive anount of 
experimental error then the initial guess will not be close 
enough. This problea is overcome by teiporarily sinplifying the 
■odel. 

The aodel is sirplified by elininating the foundation. The 
daaping effect that the tissue has on the bone is accounted for 
by a higher than nornal dunping in the beaa. The siiplified 
nodel has only four parameters to be determined by the SIDA. 
They are the bending stiffnes El, the fundasental frequency oj , 
and the daiping ratio of the beae and the stiffness k, of the 
spring. This version of the lodel will be referred to as the 
f our-naraaeter aodel (hPH). 

A reduction in the nuaber of parameters in the lodel is 
accompanied by a reduction in the nuaber of ainiiun and laxinui 
points in the error function. This increases the chance for 
convergence to the correct solution when applying the SIDA. The 
results fros applying the SIDA to the 4PH are used as part of an 
improved initial guess for the 6PH, thus increasing the chance 
for convergence when applying the SIDA to the 6PB. The process 
described herein occurs iu three phases. The SIDA is applied in 
a different way in each phase. 

In phase odo, the SIDA is applied to the 4PH. The initial 
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guess is determined by solving equations (4.7), (4.8) , (4.28) 


and (4.29) for the four parameters 

El = a*bV3L$ Zu*vPu~ (P*o/Pm,„> 2 / [(P^x/Pm.n) 2 - 1] (5.9) . 

« ■ Pmin (5.10) 

^ ZmhiPmdj/Z^P^ { (P^/P MlC ) 2 ” 1 ] / (Pm4|/Pmim) 2 (5. 11) 

k * 2^>wP uwu (Pn^/Pm, w ) 2 (5.12) 


where (Pusw f » (Pmho^Z^u) and ( Pmaa » — mj») are the key data points 
as shown in Figure 4.20, and $ is a constant which depends on 
the boundary conditions of the beam (see Table 4.2). 

In phase two, the SIDA is applied to the 6PH. However, only 
the foundation parameters are allowed to vary. The bending 
stiffness and the fundamental frequency of the beam and the 
stiffness of the spring are held constant at the values 
determined from phase one. The damping ratio of the beam is 
reduced to the reasonable value of 0.05 as mentioned earlier. 
Phase two allows the values of the foundation parameters to be 
improved without disturbing the beam and spring parameters. The 
initial guess is partially based on experience with simulating 
DPHI data "by hand" and partially based on equation (4.16). The 
fundamental frequency and damping ratio are guessed from 
experience to be one-half of those >f the beam. The mass per 
unit length is calculated by solving equation (4.16). Thus the 
initial guess is calculated from 


p ; «= ir*/L* 2(£ - 0. 05) /I (]f/2) n 


I 


iOf * uj/2 

Tr - r/2 


(5.13) 

(5.14) 

(5.15) 


O 


In phase three, the SIDA is again applied to the 6PH. All 
six paraaet« rs are allowed to vary. The initial guess is simply 
the results of phases one and two. The SIDA converges to the 
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correct solution for the DPMI data fros almost any reasonable 
forearo or leg vibration test. Examples mill be given in 
Chapters 71 and VII. 

THE COMPUTER PROGRAM 


1 . The Program 

A Fortran computer prograa was written to carry out the 
process described in the last section. Due to the complexity of 
the DPMI functions being evaluated, the program is written in 
double precision. A listing of the program is given in 
Appendix F. 

The computer program is divided into three phases of the 
total process. Each phase is similar in structure. A general 
flow chart of the program is shown in Figure 5.2 and a more 
detailed flow chart of one phase is shown in Figure 5.3. Control 
passes through the tain loop of each phase of the program until 
the iterations are terminated by the passing of one of the four 
tests as indicated in the diamond shaped boxes in the flow 
chart. 

The first test is to determine whether or not a negative 
value was obtained for one of the parameters in the previous 
iteration. Onlike the other three tests, the consequence of 
passing this test depends on the phase. In phase one, the 
parameters are returned to their old values. In phase two, the 
tissue parametric values are returned to their initial guess. In 
phase three, the 6PM is disregarded and the parametric values 
obtained for the 0PM are recalled. 

- 11*1 ii ■ I I H dii i '■ -- - - 
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The second test is to deternine whether or not the value of 
the error function has increased in the last iteration. If it 
has, then this is an indication that the paraaetric values are 
either loving away from the correct ciniaui point of the error 
function toward a aaxiaua point or that the SIDA has over- 
stepped the minima point. In either case, the old set of 
paraaetric values are closer to the correct solution than the 
new set. Therefore, the paraaeters are returned to their old 
values. 

The third test is to deternine whether or not convergence 
has occurred. Convergence is considered to have occurred when 
all of the percentage changes in the paraaeters have becone less 
than one-tenth of a percent. 

The fourth test is to deternine whether or not ten 
iterations have occurred. A liait of ten iterations is placed on 
each phase to insure that the iterations will not go on 
indefinitely. 

If all four tests fail in a given iteration, then control 
is transferred back to the top of the loop and another iteration 
is carried out. 

2. The flatrix Equation 

Within each iteration of the SIDA, a latrix equation of the 

fori 

[A] [APJ * (B) (5.5) 
is generated. The solution is to be obtained within the computer 
progaa using the subroutine DGELG froi the IBH Scientific 
Subroutine Package (*SSP) . DGELG solves the latrix equation 
nsing Gaussian eliaination. 
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Accuracy of the calculations is an important factor since 
it can influence convergence of tbe SIDA. Hatrix A, however, is 
an ili-foraed matrix, i.e., its elenents vary in absolute value 
as much as ten to twenty orders of magnitude. Ill-forned 
aatrices are very difficult to solve accurately. Therefore, 
eguation (5. 5) will be aodified to eliminate tbe ill-fornedness 
of natrix A. 

Consider natrix eguation (5.5) in conponent forn 

M 

2 Ay APi » Bj ; j * 1,2.. .H (5.16) 

fl 

Eguation (5.16) represents H linear algegraic eguations, where H 
is the number of parameters to be deternined by the SIDA. Each 
of the algebraic eguations can be nultiplied by a constant 
without altering the solution. 

M 

^ (C, A ) AP i * (CjBj) ; j = 1,2. ..H (5.17) 

i*l 

where Cj , j = 1 , 2 ...H is a set of H constants. Furthernore, the 
coefficients of each unknown can be nultiplied by a constant if 
that unknown is divided by the same constant. Using the same set 
of H constants, the symetry of matrix A is preserved. 

M 

^ (C t C, Ai j ) (&P t /C t ) = (C;Bj) ; j= 1 ,2. . • H (5.18) 

thus the new matrix eguation is 

[A] {£>} = [B] (5.19) 


where 

A t j = Cj C j A^j 

£P t = AP t /Ci (5.20) 

Bj « CjBj 

Befering to the definitions of A t j and Bj given in 
eguations (5.6) and (5.7), the orders of magnitude of each of 
the quantities in eguation (5.16) are as follows: 

A t j has order P“* P 7 1 


AP t has order P t 
Bj has order Pj*. 

The difference in the orders of magnitude of AP^, i * 1,2... tl 

results in the ill-f ornedness of natrix A. However, natrix £ 
will be well-foned if the constants are chosen so that each 
element of aatrix A is of order one. This can be accomplished by 
choosing 

C, = 1/B t ; i * 1,2.. .H (5.21) 

Then the new aatrix eguation becomes 

[X] (Ap) = (I) (5.22) 

where 

A t j = A v j /( B t B j ) AP t = B t AP L (5.23) 

and all of the components of the column aatrix I are unity. 

Matrix equation (5.22) can be solved without loss of 
accuracy because aatrix i is well-foraed. The solution, however, 
is different from the solution to aatrix eguation (5.5). The 
relationship between the two solutions is known from equation 
(5.23). Hence the solution to eguation (5.5) is calculated fron 
the solution to eguation (5.22) by 

AP L = AP t /B t ; i = 1,2. ..fl (5.24) 

3. Input 

To lake the computer program user oriented, the input 
required to run it has been simplified as Buch as possible. Only 
four lines of inforaation are required in addition to the data 
points themselves. The input is checked by the coaputer program 
and error messages are printed out to inform the nser if it is 
not in proper form. An example of input is given in Figure 5.4. 

The first line is a title. The user can insert anything he 
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wishes with a limit of sixty characters. The title is printed on 
both the output and the plot. 

The second line contains the support length of the forearm 
and the length-to- probe location ratio. This ratio is the 
distance between the left support and the driving point divided 
by the support length. The ratio must be a number between 2ero 
and one. If it is not, then an error message is printed. The 
length and ratio are read in free format. 4 

The first two columns of the third line contain an 
integer. 5 A negative integer indicates that the specimen is an 
ulna and a positive integer indicates that the specimen is a 
tibia. Recall that the boundary condition on the foundation of 
the model is either free or fixed depending on the type of 
specimen being represented. This is the only indication given to 
the program concerning the type of specimen. The data is 
interpreted according to the value given on this line regardless 
of what information is entered in the title. If a zero appears 
on this line, then the foundation is not included in the model. 

The fourth line contains the number of the data points. 
This number must also appear as an integer in the first two 
columns. It least eight but no more than sixty data points are 
allowed. An error message is printed if this is violated. 

Starting with line five, the remaining lines contain the 
data points, one per line. The forcing fregaency, magnitude of 


4 Free format: There is no restrictions on the form of the 

number, i.e., with or without a decimal point, with or without 
scientific notation. A comma and/or at least one space must 
appear between each entry. 

* Integer: Decimal points and scientific notation are not 

allowed. Rote, a one -digit number with no sign must appear in 
column two. 
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the DPHI and the phase angle of the DPMI Bust appear in order 
and in free foraat. 

The only other restriction on the input concerns units. 
Frequencies and phase angles are entered in Hertz (cycles per 
second) and degrees, respectively. All other quantities Bust 
have consistent units. Ho conversion factors have been written 
into the prograa. The CGS systea is suggested, i.e., all 
quantities are expressed in teras of centimeters, graos, seconds 
and dynes. 

JL. Output 

To take the prograa user oriented, the output aust be easy 
to read and interpret. An exanple of output is given in Figure 
5.5. The corresponding coaputer plot is given in Figure 5.6. 

The title, given by the user on the first liDe of the 
input, is printed at the top of the output page followed by the 
length aDd ratio. The paraaeters of the nodel are listed with 
their values. The data points and their corresponding DPMI's of 
the aodel are tabulated. Finally the value of the error function 
is given. 

A cocputer plot is also generated as part of the output. 
The squares represent the DPMI data points. The solid line 
represents the DPMI of the aodel, calculated using the final 
paranetric values, deterained by the SIDA. Both the aagnitode 
and the phase angle of the DPMI are plotted to visualize the 
quality of the simulation. 


CHAPTER VI 


V ERIFICATION OP THE MATHEMATICAL MODEL 


A. IN VITRO MONKEY EXPERIMENTS 


1 . Pr opose d E xperiments 

A series of experiaents was proposed by Orne and Mandke 
(1975) to verify the ■athenatical model. These experiaents are 
designed to isolate the effects of the various components of the 
vibrating foreara systea. The experiaents involve the 
application of the test procedure, described in Section I.E., to 
a monkey ara under three different conditions. 

Th -"atony of the an and foreara of a lonkey is quite 
siailar that of a huaan ara and foreara. There are, of 

course, soae ninor differences but the similarity is strong 
enough so that the results of these experiaents will provide and 

indication of the validity of the application of the 

/ ■ 

aatheaatical aodel to experiaents done with either spjcies. 

A fev aodif ications, including the addition of a fourth 
condition, were introduced before the experiaents were conducted 
by Peterson (1977). A description of the experiaents (in 
aodified fora) is given here. 

The ara of a sacrificed aonkey is disarticulated at the 
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shoulder and immediately frozen to naintain freshness until the 
experinents could be perforaed. The specinen was thawed and 
allowed to cooe to rooa temperature before testing. The 
following experiaents were then perforaed as quickly as 
possible. 

The aonkey are is positioned in the test fixture. A weight 
is placed at the top of the huaerus to represent the dovnvard 
force applied through the huaerus by the live subject, as shown 
in Pigure 6.1. This first condition should resemble an in vivo 
test as auch as possible. The driving-point aechanical iapedance 
(DPHI) of this systen is aeasured. 

A small piece of skin is reaoved from the forearm to allow 
the probe to be applied directly to the ulna. This is the second 
condition. The DPHI is again aeasured. 

All of the tissue surrounding the bones between the 
supports is reaoved. The joints and the tissue surrounding the 
joints at the supports is left intact. Care is taken that the 
support conditions are not altered between the first three 
conditions. A third set of DPHI data is taken. 

Finally, the ulna is completely excised. Holes are drilled 
in the ends of the ulna to accommodate small steel pins. Care is 
taken in drilling the boles so that the orientation of the ulna 
is not changed between the third and fourth conditions. The pins 
are supported in brackets as shown in Figure 6.1. The fourth set 
of DPHI data is taken. 

2i The Hath eaatical Hodel 

The DPHI plots produced by the experiaents are to be used 
to verify the aatheaatical aodel described in Section II. E. To 
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do this, the DPMI plots are siaulated using the aatheciatical 
aodel in its appropriate fora. The validity of the aatheaatical 
■odel is verified by deaonstrating its capability to accurately 
sinulate each of the DPHI plots produced by the ezperiaents. 
Furtheraore, each paraoetric value obtained by the sinulations 
Bust be within a range of reasonable values and, of course, aust 
be non-negative. 

In the fourth condition (excised ulna) , the ulna is 
supported by a pin and bracket at each end. The pins, which are 
Bade of steel are saootb and relatively rigid. The sooothness of 
the pins produces essentially no resistance to rotation while 
their rigidty provides essentially infinite resistance to 
translation.* Therefore the excised ulna can be aodeled as a 
sinply-supported beau. For each successive condition, in 
reversed order, the element is added to the nathenatical aodel 
which corresponds to the component of the systea which was 
reaoved in obtaining the previous condition. 

The third condition (ausculature reaoved) differs froa the 
excised-ulna condition only in the Banner in which the ulna is 
supported. Ideally, the joints provide siaple supports for the 
ends of the ulna, yielding identical DPHI plots for the two 
conditions. If the two DPHI plots are not identical, however, 
then the DPHI plot cf the an in the ausculataro-renoved 
condition will provide and indication of the true boundary 
conditions of a live foreara. 

The second condition (probe on ulna) bas all of the tissue 


• The relative rigidity of the pins vas verified by calculating 
the static stiffness of a pin and comparing it to a typical 
value of static stiffness of a bone. A difference of two to 
three orders of aagnitude vas foond. 
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surrounding the ulna and radius in place. The layer of skin 
between the probe and ulna in this case has been reaoved. 
Therefore the aatheaatical aodel includes the foundation but not 
the spring-in-series. 

Finally, the first condition (intact ara) is aodeled 
according to the aatheaatical aodel description given in Section 
II. C. Since all of the coaponents of the vibrating foreara 
systen are present, all of the eleaents of the aatheaatical 
aodel are present. 

The fora of the aatheaatical aodel for each successive 
condition (in reversed order) contains all of the paraaeters 
present in the previous condition together with one or aore 
additional paraaenters. The paraaetric valnes obtained for the 
previous condition are preserved while values for the additional 
paraaeters are obtained using the systeas identification 
algorithm (SID1) described in Chapter V. This consistent 
building-block approach to aodeling the intact ara gives greater 
confidence that the aodel actually represents e physical 
systea and that arbitrary curve-fitt; ng is reduced to a ainiaua. 

3 . Appl ic at io n of the Systeas Ident .i f ication Alqoritha 

A set of coaputer prograas was written to carry out the 
siaulations discussed above using the SIDA. These coaputer 
prograss are each si, : lar to a "one-phase" version of the 
coaputer prograa described in Section T.D. The aost significant 
aodification is that the derivatives, calculated within each 
iteration of the SIDA, are replaced by finite differences, 

X • © • 0 


dZn/dP, - [ Z n ( P*, ♦ tP, ) - Z n (PJ]/AP 


( 6 . 1 ) 
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where the finite increment in the parameter 6P t , is taken to be 
one percent of the current value of the parameter. 

Early on in the development of the SIDA, some sets of data 
were simulated using the SIDA both with exact derivatives and 
with the derivatives approximated by finite differences. The 
values of the derivatives and the finite differences were 
printed out by the computer programs so that they could be 
compared. Their values were found to be in agreement within at 
least two, and often within three decimal places. Hence, 
accuracy of the finite differences does not present a problem. 

A trade-off exists between the effort spent in deriving 
exact expressions for the derivatives of the DPKI function with 
respect to each parameter of the mathematical model and computer 
time spent in calculating the finite difference approximations 
to those derivatives. The set of computer programs used to 
simulate the in vitro experiments must be very adaptable. 
Several versions of t he mathematical model are use in an attempt 
to produce good simulations, but each corresponding version of 
the program is run only a few times. Hben finite differences are 
used rather than exact derivatives, much less effort is reguired 
to change the computer program and employ a different version of 
the mathematical model. Therefore, the extra computer time spent 
to calculate the finite differences is justified by their 
adaptability aud convenience. On the other hand, the coapater 
program which was developed to simulate in vivo tests is to be 
run many times without changes. The same program is used with 
many different sets of data* The effort spent in deriving exact 
expressions for the derivatives reguired for this computer 
program is justified by the saving of much computer time. 
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Another important difference is that the set of "one-phase" 
co a pa ter programs is not as user-oriented as the coapater 
program described in Section V.D. Adaptability is required not 
only in the mathematical model bat also in the method of 
establishing the initial guess. Therefore, the initial guess is 
calculated "by hand M and read in at the beginning of the 
computer program. This adaptability is more important than the 
simplicity of tha input in this case. 

Q. Results From Monkey 663 

The series of experiments, described earlier in this 
section, were performed on the forearms of three monkeys, 
identified by their numbers, 659, 663 and 665. The DPMI data 
produced by these experiments were simulated by the set of 
computer programs discussed above, using the various versions of 
the mathematical model. The resulting DPMI plots associated with 
Monkey 663 are shown in Figures 6.2 through 6.5. The solid lines 
represent the DPMI of the mathematical model while the boxes 
represent the data points generated in the experiments. The 
corresponding parametric values are listed in Table 6.1. 

Figure 6.2 is the DPMI plot of the nlna in its excised 
state. As expected, the DPMI data is veil simulated as a simply- 
supported beam. Therefore the value obtained for the bending 
stiffnesses the best possible estimate of its true value. 

Figure 6.3 is the DPMI plot of the same ulna with the 
musculature removed but with the joints left intact. It is 
easily seen that the excised-ulna and musculature-removed plots 
are quite different from one another other. The excised ulna is 
virtually simply-supported. Therefore, since no other parameters 

c - 
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were changed between the excised-ulna and musculature-removed 
cases, the support conditions of the ulna when the joints are 
intact Bust be something other than simply-supported. 

The bone parametric values determined fron the excised-ulna 
case were used for the musculature-removed case; holding then 
constant while determining values for the boundary condition 
paraneters that will best siaulate the data. The boundary 
conditions which produce the best results were found to be a 
rotational spring on one end of the bean and simply-supported on 
the other. Damping was also included at both ends of the beam. 

Based on the parametric studies of Section IV. B, a 
significant amount of resistance to rotation can be created if 
the downward force applied through the humerus is not directly 
in line with the support as shown in Pigure 6.6. Hence, this is 
most likely the major cause of the resistance to rotation at the 
support in these experiments, but experimental verification is 
necessary. 

Figure 6.4 is the DPHI plot of the arm in which the layer 
of skin between the probe and ulna is removed but the rest of 
the tissue is left intact. The major difference between the 
musculature-removed and probe-on-ulna plots is the increase in 
damping in the latter case, i.e. , the region around the minimum 
point of the DPHI plot is moved upward. The tissue, in fact, 
contributes much more to overall damping than does the bone. 

Pigure 6.5 is the DPHI plot of the intact arm. The major 
difference between the probe-on-ulna and intact-arm plots is an 
overall decrease in DPHI. This is to be expected since the skin 
between the probe and the bone is in series with the bone. The 
DPHI of the whole system is less than the DPHI of either part 
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alone. 

A slightly better fit is obtained using Orne’s three- 
parameter aodel for the shin (see Section II. C) rather than the 
spring alone. However, the skin, when tested alone, does behave 
as a sinple spring, see Figure 2.2. These experiments would have 
to be rerun to include higher frequencies to better define this 
behavior. 

Since the mathematical aodel has all of the capabilities 
necessary to simulate the entire set of in vitro experiments, it 
is a good representation of the physical system. In dealing with 
an in vivo test, however, the support conditions of the physical 
system must be evaluated. The parametric values obtained from a 
simulation in this case, will be valid only if the boundary 
conditions of the aathematical aodel are a good representation 
of the support conditions of the physical system. 

5. Pesults From Bonkey 665 

The experiaents run on Bonkey 663, as discussed above, were 
also run on Bonkey 665. The data was simulated using the SIDA 
and the sane versions of the aathematical model. The resulting 
DPMI plots are shown in Figures 6.7 through 6.10. The 
corresponding parametric values are listed in Table 6.2. 

Again, the excised-ulna data of Figure 6.7 is well 
simulated as a simply-supported beam. The remainder of the data 
sets, however, are not siaulated as well. A disturbance, 
occurring at about 200 fir, in the ausculatnre reaoved plot 
becomes progressively more pronounced in the probe-on-ulna and 
intact-arm plots. This disturbance is similar in appearance to 
that which is expected fron the tissue surrounding the bone. 
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However, the tissue is not the cause of the disturbance in this 
case since it also appears in the musculature-removed plot. The 
true origin of the disturbance in this data set is not known. It 
is not likely, however, that it is a true characteristic of the 
vibrating forearn system, since it does not appear in the data 
froa the other two monkeys. 

The disturbance found to occur in nost of the data from 
human sabjects is still thought to be a result of the tissue 
surrounding the bone. This situation does not occur in the 
monkey data, since the monkey has less tissue on his bones. 
Similar experiments on a human cadaver arm must be run to verify 
this effect. 

6. Results From Monkey 659 

The DPMI plots of Monkey 659 are shown in Figures 6.11 
through 6.15. The corresponding parametric values are listed in 
Table 6.3. Two major differences exist between the procedure of 
these experiments and that of Monkeys 663 and 665. First, DPMI 
data for the ulna in its excised state was not taken until two 
months after the other DPMI data. During that time, the ulna was 
stored in a refrigerator. Second, DPMI measurements were taken 
on the intact arm at both a 400 and 600 gram-force preload. 

As in the other two cases, the excised-ulna data of Figure 
6.11 is , well simulated as a simply-supported bean. The 
musculature-removed data of Figure 6.12, however, could not be 
simulated directly using the same boundary conditions in the 
mathematical model as those used for Monkeys 663 and 665. Becall 
that the excised-ulna data was obtained two months after the 
other data. Although the attempt was made to maintain freshness. 


significant deterioration had occurred. In fact, the SIDA 
indicates a thirty-two percent decrease in the bending stiffness 
of the ulna over that tine, with this change in bending 
stiffness taken into account, a good simulation was obtained for 
the musculature-removed plot. 

The probe-on-ulna data of Figure 6.13 is well siaulated by 
the mathematical model. 

Figures 6.14 and 6.15 are DPMI plots of the intact ara with 
400 and 600 gran-force preloads, respectively. As with the data 
from Monkeys 663 and 665, the presence of the skin between the 
probe and bone has the effect of decreasing the DPMI. This 
decrease is less for the 600 than for the 400 gran-force 
preload, as expected. If the preload could be nade high enough 
without destroying the ulna, the decrease in DPMI would 
eventually disappear altogether. 

Jk BEHDISG TESTS 


1 ♦ Procedure 

The DPMI technique and its analysis described herein, 
results in a value measured for the bending stiffness of a long 
bone. To verify that this neasurenent is valid, the bending 
stiffness of an excised long bone, which has been aeasured using 
the DPMI technique, was measure d using another independent 
technique. Each technique should give the sane result. The 
alternate technique involves a sinple three-point bending test 
fron which a load-deflection curve is generated. 

The ulna of Monkey 659 was tested in each of four 
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conditions described in the last section. After the tests were 
completed, it vas wrapped in gau 2 e, soaked in Einger*s solution 
and frozen to maintain as much freshness as possible. The ulna 
was then nailed from Stanford University, California to Vayne 
State University, Hichigan, where it was again frozen. Just 
prior to testing, the ulna was brought to room temperature by 
soaking it in a jar of Binger*s solution. 

h Haterial Testing System (UTS) machine was used to perform 
the bending tests. The ulna, already pinned from the !>PHI test, 
was placed in the bending fixture as shown in Figure 6.16. The 
UTS machine was programmed to apply a constant deflection rate 
to the center of the ulna. Several different deflection rates, 
ranging from 0.5 x 10 -3 to 0.5 in/s (1.27 x 10” 3 to 1.27 cm/s) 
were used. These deflection rates are slow enough so that mass 
and damping effects are not present. The maximum deflection, 
approximately one-half centimeter, produced stresses which are 
within the elastic range. The load-deflection curve, shown in 
Figure 6. 17, was generated on an x-y recorder, using the force 
and displacement signals from the ATS machine. 

The static stiffness of the ulna is determined from the 
load-deflection curve using the relation 

K = AF/A5 (6-2) 

where £F and && are shown in Figure 6.17. The bending stiffness 
of the ulna, using a uniform beam model, is then determined from 
the relation 

El * KL*/48 (6.3) 

The ulna was allowed to dry for a period of two months. The 
value of the bending stiffness was then remeasured. 
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2. Resalts and Evaluation 

A sueaary of the Measurements described above is given in 
Table 6.4. Note that the value obtained for the bending 
stiffness from each successive test is significantly lover than 
that obtained froa the previous test. 

Although the atteapt vas aade to keep the ulna as fresh as 
possible, it had deteriorated to soae degree. Table 6.4 suggests 
a trend towards lover values of bending stiffness as the ulna 
deteriorates. Therefore, higher values vould be expected if the 
bending tests had been perforaed iaaediately after the DPHI 
tests. The percent difference vould then be reduced, if not 
eliainated all together. 

With the effect of deterioration taken into account, the 
bending stiffness values aeasured by the tvo independent 
techniques are fairly consistent. Hence, these results support 
the validity of the DPHI tests. 

£. NON-BIOLOGICAL t ests 
1 . The Sy st eas 

To verify that the eguipaent is actually aeasuring the DPHI 
properly, the DPHI of tvo non-biological systeas is Measured. 
Non-biological systeas can be constructed in such a vay that 
their aechanical response is auch aore predictable than that of 
a biological systea. Furtheraore, the coaponents of that systea 
can be aade of aaterials vhose aechanical properties are veil 
known. In particular, the tvo systeas at hand are aade of coaaon 
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The first system is siaply the calibration aass discussed 
in Section I.E. The second systea consists of a unifora bean, 
aachined fron a bar of aluninua, and supported by pins near its 
ends. 

2. Calibration Bass 

The calibration aass is cylindrical in shape, is Bade of 
brass and has a aass of 98.4 grans. The aagnitude of the DPBI of 
a pure aass is (see Table 2.1) 

Z = ap (6.4) 

Therefore, a log-log plot of the DPBI data should fora a 
straight line on a +45° angle. However, this is true only for 
relatively low frequencies. It very high frequencies the aass 
deforos. Therefore, the DPBI curve should go through a series of 
resonant and anti -resonant points. 

DPBI data, taken for the calibration aass up to a frequency 
of 3000 Hz, is shown in Figure 6.18. The calibration aass 
?ibrates as a pure aass up to a frequency of about 1000 Hz. It 
then approaches its first anti-resonant point at approximately 
2800 Hz. 

The systea is aodeled as a simple one-diaensional 

continuous rod with a haraonic force applied to its base. The 

* 

DPBI of such a aodel is 

Z* = inp tanV> / (6.5) 

where 

r * Tip/to 

a* aass 


<*> * fundamental frequency 


p = forcing frequency 

The aathematical model was used with the SIDA, described in 
Chapter V, to determine that the fundamental frequency of the 
system is 5519 Hz (i.e., the first anti-resonant point is u/2 * 
2760 Hz) . The DPHI of the model is shown as a solid line in 
Figure 6.18. 

The modulus of elasticity and density of brass are known 
quantities aud the height and diameter of the calibration mass 
are easily measured. The fundamental frequency , estimated from 
to = tr/L fiTF (6-6) 

is found to be on the order of 60,000 Hz, with a corresponding 
anti-resonant point at 30,000 Bz. 

Since the anti-resonant frequency determined from the DPHI 
data is a whole order of magnitude lower than the expected 
value, it must be a sub-anti-resonant 7 point. If the frequency 
range of the DPHI data could be extended beyond 3000 Hz, sub- 
resonant and more sub-anti-resonant points would be observed. 
These points may be due to the deformation of the screw 
connection between the impedance head and the calibration aass. 

3. Aluminum Beam 

The aluminum beam and its support brackets are shown with 
their dimensions in Figure 6.19. The purpose of the aluminum 
bean is to provide a standard to insure that the impedance 
equipment is operating properly each tine it is used. Hetal, 
unlike biological materials, remains unchanged over a long 


v "Sub-anti-resonant" point refers to a local disturbance whose 
source is outside the system of interest; analogous to "sub- 
resonant" point, see Section IV. D. 
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period of tine. Therefore the true DPMI of the aluninun beaa 
vill remain unchanged. The DPMI plot of the aluninun beam should 
be generated prior to each use of impedance eguipnent. If any 
deviation appears in this data then the eguipnent should be 
checked for aalf auctioning. 

The aluninun bean vas designed to have a static stiffness 
and fundanental frequency in the sane range as a typical aonkey 
ulna. Unfortunately, it is not possible to produce a unifora 
bean with these properties and with a cross section large enough 
to accoaodate rigid support pins. Therefore, it vas necessary to 
sake the ends of the beaa larger in cross section than the aid- 
portion. Only a very snail effect on the DPMI data plot due to 
the enlarged ends is anticipated. 

A typical set of DPHI data fron the aluninun beaa vas 
simulated using unifora, sinply-supported beaa aodel vith the 
SIDA. Its DPMI plot is shovn in Figure 6.20. 

The aodulus of elasticity E, and density p, of aluainua is 
known 8 and the dimensions of the beaa are given in Figure 6.19. 
The bending stiffness and fundanental frequency are calculated 
using 

El « Evd*/64 (6.7) 

to ■ (n*/L*) jEI/pA * (TT*d/4L*) {I/p (6.8) 

An area-aoaent aethod of analysis and a Bayleigh aetbod analysis 
vere carried out to determine the effect of the enlarged ends on 
the bending stiffness and fundanental frequency, respectively. 
These values, together vith those determined froa the DPHI data 
are listed in Table 6.5. 


• E * 7x10»» dyne/ca*, p * 2.7 g/ca*, e.g., see Faires (1965). 


It is seen fro* Table 6.5 that significant differences are 
apparent between the predicted values of the bendurj stiffnes 
and fundamental frequency and their values deternxned froa the 
DPHI data. So»e, but not all# of that difference is accounted 
for by including the effect of the enlarged ends of the bean. 
The only other possible source of the discrepancy (assuoing# of 
course# the inpedance equipaent is functioning properly) is in 
the boundary ' conditions. It was shown in Section IV. B# that 
resistance to rotation at an otherwise siaple support of a bean 
tends to aove the DPHI curve upward and to the right. Therefore# 
there night be and excessive aaount of friction in the pins 
which support the bean. 1 light oil should be applied to the 
pins to eliainate this friction. 

It is seen fro* Figure 6.20 that the DPHI data and the DPHI 
of the *athe*atical aodel are well correlated up to a frequency 
of about 1000 Hz. The anti-resonant point of the aathenutical 
■◦del# however# is a few hundred Hz higher than the anti- 
resonant poiu . of the systen. 

Recall that a sob-anti-resonant point vas observed near 
this frequency in the DPHI data of the calibration nass# aost 
likely due to deforaations in the screw connection at the 
inpedance head. It is possible that a siailar sob- anti-resonant 
point is occurring due to the screw connection between the 
inpedance < bead and the probe. This sub-anti-resonant point say 
or aay not be exactly the saae frequency as the previous one. 
Since the sub-anti-resonant point is relatively close to the 
anti-resonant point of the beaa# the observed anti-resonant 
point is a combination of the two. 
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CHAPX£B TXX 

APPLICATION TO EXISTI NG DATA 


it THOMPSON 1 S ORIGINAL DATA 


1 . Results 

Thoapson tea so red the driving-point aechanical iapedance 
(DPHI) of the foreara of several hoaao subjects using the 
icpedance aeasuring eguipaent which he developed (see Section 
I.E). The tests were perforaed over a freguency range froa 50 to 
1000 Hz using three different preload forces. The systeas 
identification algoritha* (SIDA) was then used to deteraine 
paraaetric values for the aatheaatical aodel which best siaclate 
the data for eight of these subjects. The DPHI plots froa one of 
these subjects. Subject TT, are shown in Figures 7.1, 7.2 and 
7.3. The solid lines represent the DPHI of the aatheaatical 
aodel while the boxes represent the data points generated by 
Thoapson. Vh< DPHI plots of the reaaining seven subjects are 


• The coaputer prograa which incorporates the SIDA is siailar to 
the one presented in Section f.3. Tbe only difference is that 
the coaputer prograa used here has the capability of siaulating 
three sets of data siaultaneousij , '.bus deteraining the three 
values for the spring-in-series; one corresponding to each 
preload. 
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presented in lppendix 6. The corresponding paraaetric values for 
all eight subjects are listed together with other available 
information in Table 7.1. 

2. Discussion 

The paraaetric values listed in Table 7.1 are good 
approximations of the geoaetrical and aaterial properties of th? 
physical system provided the aatheiatical model is a good 
representation of that physical systea. Therefore, to 
investigate the validity of these values, it is necessary to 
evaluate the support conditions. Ill other aspects of the 
aatheiatical nodel were shown in previous chapters to oe very 
good approxiaations of the vibrating forearn systea. 

When positioning the subject’s forearm in the fixture, 
Thompson carefully lined up the huaerus with the support "by 
eye". The misalignment (discussed in Section 71. A) may not be 
perfectly eliminated but it is certainly significantly reduced. 
Thus the supports are relatively free froa rotational 
resistance. 

Thompson Bade the supports as rigid as possible mith 
respect to translation by forning plaster pads under both the 
wrist and elbow. He demonstrated the rigidity of the supports by 
shoving that the DPHI was independent of both the clamping force 
at the wrist and the downward force applied through the humerus 
et the elbow. 

Based on the discussion above, the supports are virtually 
simply-supported. The paraaetric values listed in Table 7.1 were 
obtained using the S1D1 and the aatheaatical model with simply- 
supported boundary conditions. Therefore these values are very 
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good approximations to the actual geometrical and material 
properties of the vibrating forearm system for each subject. 
Furthermore, the simulations appear to give accurate results. 
The error functions associated with each plot are within two 
percent and about half of them are within one percent. 

The mass per unit length of the bone p, is calculated by 
solving equation (3.4) 

p = (tt/L)* EI/a* (7.1) 
where L is the support length, El is the bending stiffness and to 
is the fundamental frequency. This value represents the total 
mass per unit length of the bone. Measurements of bone mineral 
content (BBC) were also taken for each subject using a Norland- 
Cameron Bone Hineral Analyzer. This value represents the mineral 
mass per unit length of the bone. Values for p and BBC for each 
subject are alsc listed in Table 7.1. It is reasonable to expect 
these two quantities to correlate quite well since all bones 
tested are bones of healthy, young adults. The correlation 
coefficient r is in fact 0.81, a reasonably high value. 

Strong correlations have been found to exist between 
bending stiffness ana BBC. (see Borders, Peterson and Orne, 1977 
and Jurist and Foltz, 1977). Since the existence of this 
correlation is well established, it is reasonable to expect a 
similar correlation between the values of bending stiffness and 
BBC listed in Table 7.1 provided the values for bending 
stiffness are valid. The correlation coefficient r, of such a 
correlation, was found to be 0.87, a value comparable to 
findings of Borders, Petersen and Orne (1977) and Jurist and 
Foltz (1977). 

Each of the points discussed above support the validity of 
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parametric values listed in Table 7.1. 
B. Ionket data 


1 . Results 

Since the developaent of Thompson's impedance Measuring 
equipment, DPHI data has been generated on a routine basis for 
the forearas and legs of aonkeys at Aaes Research Center. 
Ninety-four sets of such data froa twenty-six different monkeys 
have been Bade available through personal communication. These 
tests were run over a frequency range froa 100 to 2000 or 3000 
Hz. The preload force in aost cases was 600 graa-force (589 x 
10 3 dyne) , although soae tests were run with both a 600 and a 
300 gran-force (294 x 10 3 dyne) preload. 

The computer program presented in Section ?.D was used to 
determine parametric values and generate a DPHI plot for each of 
these sets of data. A representative set of six of these DPHI 
plots are presented in Figures 7.4 through 7.9. They are froa 
the tests run on the leg and foreara of Bonkeys 2, 16 and 17. 
The corresponding parametric values are listed with other 
available inforaation in Table 7.2. 

« 

2. Discussion 

The DPHI plots of the legs appear to indicate that the 
siaulations are guite accurate. The DPHI plots of the forearas, 
however, indicate that the siaulations are not accurate. 
Furtheraore, the SIDA did not converge when applied to two of 
these data sets (forearas of Bonkeys 2 and 16) using the six- 
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parameter model* 0 (6PB) . This trend is present throughout most 
of the data. 

To investigate the validity of the parametric values listed 
in Table 7.2, it is necessary to perforn two evaluations. First, 
the cause of the difference between the leg data and the foreara 
data must be deterained. Second, the support conditions aust be 
evaluated. 

Fora the ratio k/K using the paraaetric values listed in 
Table 7.2, where 

K = 3EIL/a*b* « the spring constant of the bone 
= 48EI/L* for the probe at the center (tibia) 

= 625EI/12L 3 for the probe at ,6L (ulna) 
k * spring constant of the skin 

The value of k/K is also listed in Table 7.2 for each linb. 
Since k/K is the ratio of the spring constants of the skin and 
bone, it is a major factor in determining the magnitudes of the 
DPMI data. The stiffness of the skin k, is Bade as high as 
possible by increasing the preload force on the electromagnetic 
shaker to a tolerable limit. If it were possible to increase k 
to infinity, then the resulting DPMI plot would be that of the 
system without the skin. If k is relatively low, such that k/K 
is equal to 2 or 3, then aost of the characteristics of the 
underlying systea will be "masked" by the presence of the skin. 
Therefore > k/K aust be high enough to "expose" all of the 
characteristics of the rest of the systea. 

In view of these coaaents, examine the values of k/K listed 


10 The results obtained from applying the SID1 with the four- 
paraaeter model (4PH) are presented in these cases. For an 
explanation of the 4PH and the 6PM, see Section V.C. 
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for each limb in Table 7.2. Note that in general, k/K is such 
higher for the foreara than for the leg. In particular, note the 
extreaely high values for the forearn of Bonkeys 16 and 17. k 
is relatively constant since the preload is the sane for all 
bones. Bence, it is reasonable to expect this large variation in 
k/K because the bone stiffness K, varies significantly with the 
size of the bone. Table 7.2 shows that K is alnost an order of 
aagnitude larger for the tibiae than for the ulnae. Therefore 
the data which best exhibits the characteristics of the systeo 
(the bone, tissue and supports) are those of the forearas 
because k/K is greater. Furtheraore, the aost revealing foreara 
data is from Bonkeys 16 and 17. 

Referring to Figures 7.8 and 7.9, it can be seen that there 
is an additional relative aininum in the DPBI data at about 1200 
Hz which the siaply-supported beaa aodel can not account for. 
This is typical of the sets of data which have a high k/K value. 
Based on the above discussion regarding the Masking effect of a 
low k/K value, it is reasonable to suspect that this additional 
relative ainiaua is characteristic of aost of the liabs but that 
it is hidden by the low k/K value in aany cases, particularly 
with the legs. 

In Section IV. B, it was shown that although the boundary 
conditions do not affect the shape of the DPBI curve at low 
freguency, they can affect it at high freguency. A DPBI curve 
with a general shape siailar to that of the DPBI data in Figures 
7.8 and 7.9 can be generated if the boundary conditions of the 
beaa are those of case 5, i.e., a translational spring (and 
daaper) at each end. This is further deaonstrated by the non-’ 
diaensionalized DPBI plot shown in Figure 7.10. This figure was 
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generated using the following non-dimensional parametric values: 

'S * 0.1 

T, = 2k,L>/EI = 10 
T z = 2k* L*/EI = 10 
C T , c c.cu/k, * 2 
C Ti * c t uj/* z * 2 
k/K = 20 

Furthermore, the masking effect of a low value of k/K is 
accounted for in this model as demonstrated by Figure 7.11, 
where its value was reduced from 20 to 2. Knowing the type of 
boundary conditions which can possibly produce the kind of DPHI 
data in Figures 7.6 and 7.9, speculations can be made on the 
cause of such data. 

it some point in the development of the impedance measuring 
procedure, the plaster pads at the supports (discussed in 
Section VII. i) were replaced by putty. Host putty exhibits both 
springlike and daoperlike behavior. Therefore, it is very likely 
that the putty is a major factor in producing the second 
relative minimum in the DPHI data. Furthermore, it is difficult 
to rigidly support the tibia at the ankle. The soft tissue 
surrounding the tibia nay also be contributing to the springlike 
and danperlike behavior of the support. 

The boundary conditions of the mathematical model used to 
obtain the parametric values listed in Table 7.2 are simply- 
supported. Since the support conditions of the forearms and legs 
for the DPHI tests discussed above are not simply-supported, the 
paraietric values are not accurate. 
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CHAPTER VIII 
CONCLOSIOH 


A. SOHHARY 


1. Overview 

A brief summary of the research project as a whole is 
given, followed by a summary of the contributions of this work. 
It is important to consider the relationship between this work 
and the work of other investigators involved in the research 
project and to give then appropriate credit. 

The iapedance measuring eguipaent and procedure were 
developed by Thompson (1973). He aeasured the driving- point 
mechanical iapedance (DPHI) in vivo of tne forearm of several 
healthy, young, adult, human subjects. Thompson also used a 
single-degree-of-f reedom oscillator (SDOFO) in series with a 
spring as a mathematical model to interpret his data. 

Orne (1974) proposed a visco-elastic beam model to better 
simulate the DPHI data. Orne and Bandke (1975) further improved 
the mathematical model to account for some of the finer details 
of the D°HI data. They also proposed a series of experiments to 
be run on a monkey forearm to verify the mathematical model. 

Petersen (1977) performed the experiments which Orne 
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proposed. ** One of the ulnae from these experiments was also 
tested statically in three-point bending on a Baterials Testing 
Systen (UTS) aachine. 

A series of experiaents involving the aeasureaent of 
breaking strength of excised canine long bones was perforned; 
see Borders, Tetersen and Orne (1977). Bending tests were 
conducted on an ATS aachine and correlations were established 
between the various paraaeters aeasnred in these tests. 

Petersen (1977) Bade soae modifications to Thonpson's test 
procedure and applied it ijn vivo to both the forearm and leg of 
monkeys. DPHI data has since been collected for monkey forearms 
aod legs on a routine basis by Howard (personal communication) 
at Ames Research Center. 

Concurrently, the mathematical model was further developed, 
in extensive paranetric study was made using the mathematical 
model. A systems identification algorithm (SIDA) was developed 
and applied to the data obtained during the experiaents and 
tests mentioned above. 

2. P arametric Study 

A parametric study has been carried out (Chapter IV) to 
determine the effect of each parameter of the mathematical model 
on its DPHI response. Two accomplishments were attained as a 
result of the study. First, an increased understanding of the 
effects of the paraaeters was gained. Second, many qualitative 
relationships between the paraaeters and the characteristics of 


»* These experiments were rerun with a wider frequency range on 
both the forearm and leg of a monkey. However, the impedance 
measuring equipment was not functioning properly and the DPHI 
data could not be interpreted. 
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the DPHI carve were derived * A brief description of the 
mathematical aodel followed by a summary of soae of the aajor 
findings is given here. 

The ulna of the vibrating foreara systen is represented by 
a uniform, linear, visco-elastic, Euler-Bernoulli beam. The shin 
and tissue compressed between the probe end bone is represented 
by a spring in series with the bean. The renaining skin and 
tissue surrounding the bone is represented by a visco-elastic 
foundation with mass. 

A linear beam aodel, regardless of its boundary conditions, 
generates a DPHI curve which is identical in shape to that of a 
SDOPO up to a frequency of at least two tines, and often as nuch 
as ten tines the fundanental frequency. This is demonstrated by 
the figures presented in Chapter IV for several different types 
of non-classical boundary conditions. The only parameter 
affecting the shape of the curve is the damping ratio. 
Furthermore, the position of the curve on the plot is entirely 
determined by the static stiffness and fundanental frequency of 
the beam. 

Hone of the boundary conditions discussed in Chapter IV 
produce a rigid body node of vibration, i.e., produce a zero 
fundamental frequency. In fact there exists only two cases of 
boundary conditions which will produce a rigid body mode: free- 
free and pinned-free. The DPHI curve in these two cases is 
identical, up to three or four times the first antiresonant 
frequency, to a SDOFO with the driving force applied to its 
base. 

A few approximate relationships between the parameters of 
the beam and the characteristics of its DPHI curve have been 
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derived. They are useful for obtaining a first approximation for 
the paraneters directly from a set of DPMI test data. 

A transfer aatrix Method of analysis vas developed to study 
the effect of taper (Appendix B) . This Method allows any 
parameter which is varying along the length of the bean to be 
approximated by a series of step functions constant within each 
eleaent of the bean. The transfer Matrix is generated from the 
exact solution of the bean within each element. (Note: the 
eguations which nake up the matrix could also be rearranged to 
form a stiffness matrix, thus producing a finite element 
representation of the beam.) 

The conclusion drawn from applying the transfer matrix 
method to a calculation of the DPHI is that the taper does not 
affect the DPHI in the frequency range of the DPHI tests. A 
uniform beam and a tapered beam with the same static stiffness 
each produce a DPHI curve which is identical op to frequencies 
of at least an order of Magnitude above the fundamental 
frequency . 

A visco-elastic foundation with mass has two effects on the 
DPHI curve of a bean. First, it produces a subresonant 
disturbance in the otherwise smooth curve. This disturbance is 
present in Many DPHI data sets. Second, the foundation produces 
a damping effect, similar to the damping in the bean. Hence, the 
minimum point of the DPHI curve is affected by the paraneters of 
the foundation. This effect could not be quantified in close 
form due to the complexity of the DPHI equations. However, 
approximate relationships were derived which are valid for sone 
range of parametric values. 

A spring in series with a bean has its Major effect in the 
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high frequency range. The total DPMI of tvo systeos in series is 
dosinated by whichever system has a lower DPMI. Thus, the spring 
dominates the total DPMI in the high fregnency range where its 
DPMI is low. The high frequency data froa a DPMI test has been 
used in the past to approxiaate the stiffness of the spring. 
However, data is not available in a high enough frequency range 
to coapletely eliminate the effect of the bean. Bence, this 
approach led tc significant errors in estiaating the spring 
stiffness, which in turn led to errors in estimating the 
stiffness of the beaa. In alternate approach has been developed 
which is such nore accurate. The approach is based on the 
location of the aaximua point of the DPMI curve which occurs due 
to the spring. This eliainates the need for the high frequency 
data, otherwise required to make the estiaate. 

3^ The Sys t ems Identification A lqoritha 

# 

A SIDA has been developed to deteraine the parametric 
values of the aatbenatical aodel which best simulate the data 
obtained froa a DPMI test (Chapter V) . The SIDA is based on 
aininizing the error function; a function siailar in fora to 
that used in a classical least-squares aethod. 

Due to the coaplexity of the DPMI equations of the 
aatheaatical aodel, the error function is very nonlinear with 
respect to its paraaeters. Consequently, a systea of equations 
obtained by setting the derivative with respect to each 
parameter equal to zero, is virtually impossible to solve. 
Bather than solving for the parametric values directly, an 
iterative procedure was developed which involves the calculation 
of a change in each paraaetric value which will bring that 


103 


paraoeter closer to its correct value. 

The expression for the DPMI of the aatheaatical aodel was 
replaced by the first two teras of its Taylor series expansion 
abont the point associated vith the current value of each 
paraoeter. Then differentiating the error function *ith respect 
to changes in the paraaetric values leads to a systea of 
eguations which are linear in these changes. To start the 
iteration procedure, an initial guess for each paraaetric value 
is obtained using the relationships derived in the paraaetric 
study. 

ft . Evaluation of Ex istinq Experiaents and Tests 

Data froo several groups of DPBI tests and experiaents have 
been aade available through personal coaaunica tion with iaes 
Research Center, haong thea are (1) in vitro aonkey experiaents, 
(2) nonbiological tests, (3) Thoopson's original in vivo tests 
and (0) aore recent in vivo aonkey tests. 

The in vitro aonkey experiaents, discussed in Section VI. 1, 
involve the aeasureaent of DPMI of a aonkey foreara in several 
stages as the ulna is being excised. The aatheaatical aodel was 
shown to be a good representation of the physical systea by 
using it in its appropriate fora to siauiate the whole set of 
experiaents vith a consistent set of paraaetric values. Bending 
tests were perforaed on one of the ulnae which were excised 
during the experiaents (Section VI. B) . These tests verify the 
value obtained for the bending stiffness of that ulna. The 
expedients, however, revealed that a problea exists in the 
consistency of the support conditions of the speciaen. This 
problea will be snaaarized In the next section. 
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DPHI tests vere run on tvo nonbiological systeas: a "rigid" 
aass and an aluoinua boa a. The data froa these tests vere 
studied. Baking use of soae simple aatheaatical aodels (Section 

VI. C). The results indicate that the iapedance aeasuring systea 
is, in fact, aeasuring the DPHI properly over aost of the 
frequency range. 

Thonpson, the developer of the iapedance aeasuring 
equipment, aeasured the DPHI in vivo of the foreara of several 
human subjects. The aatheaatical aodel vas used with the SIDA to 
deteraine the paraaetric values (Section VII. 1). The results 
indicate that both the iapedance aeasuring equipment and the 
analysis procedure are working well. Values vere obtained for 
bending stiffness of the ulna of each subject. 

The iapedance aeasuring procedure has since been aodified 
and applied to forearas and legs of aonkeys £n vivo (Section 

VII. B) . These tests revealed a forther problem vith the support 
conditions of the specimen and is also summarized in the next, 
section. 

FECOHHF HPATIOHS 

1 . Problems Revealed by Experiments 

In siaulating the in vitro experiments of Section Vl.i, 
only a fev paraaetric values vere determined froa each set of 
data. In particular, the bone paraaeters and the support 
paraaeters vere determined froa tvo different DPHI data plots. 
Bovever, when siaulating an in vivo test, values for the vhole 
set of paraaeters must be determined simultaneously froa a 
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tingle set of data. If this set of paraseters contains 
stiffnesses of both the bone and supports, then the nunber of 
paraseters will be too great. It is ispractical to use a 
sathenatical sodel which has too sany stiffness paraseters since 
it is inpossible to identify each paraseter individually. On the 
other hand, the boundary conditions of the sathenatical sodel 
sust be a good representation of the support conditions of the 
physical systes. The only way to solve this dilenaa is to have 
sooe control over the support conditions in the in vivo tests. 

Ideally, the support conditions in the in vivo tests should 
be nade sisply-supported. To do this, all sources of lateral 
translation and resistance to rotation at the supports sust be 
elisinated. k systesatic procedure should be developed which 
consistently produces support conditions which are virtually 
sisply-supported. 

In practice, it say not be possible to consistently attain 
the sisply-supported support condition. However, even if this is 
the case, a systesatic procedure is needed for positioning the 
specisen in the test fixture. Two reguiresents sust be isposed 
on this procedure. First, the procedure sust produce support 
conditions which are as nearly siaple-supported as possible (or 
pra 'tical) . The purpose in striving for such a support condition 
is to sazisize the strength of the dependence of the DPHI of the 
vibrating foreara or leg systes on the bone stiffness, thus 
saxisizing the sensitivity of the DPHI to changes in the bone 
stiffness. Secondly, the procedure sust produce support 
conditions which are repeatable. If the support conditions are 
not to be known, then they sust at least be consistent froa one 
test to the next. In this case, the value of bone stiffness 
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inferred through the aatheaatical aodel will be an index of the 
true bone stiffness rather than an absolute aeasure. 

2. Farther Suggested Experiments 

Based on the paranetric studies of Section IV. B, a 

significant aaount of resistance to rotation can be created if 

the downward force applied through the huaerus is not directly 

in line with the support as shown in Figure 6.6. It is believed, 

therefore, that this is a aajor cause of the rotational 

resistance that was found to be present at one of the supports 

in the in vitro aonkey foreara experiaents. This speculation can 

be tested by running additional in vitro aonkey foreara 

experiaents. In these experiaents, the support is to be 

positioned in several different locations in the vicinity of the 

elbow, thus varying the degree of aisalignaent . 1 value can be 

obtained for the bending stiffness of the ulna using the simply- 

supported bean aodel and the SIDi in each case. Then excising 

the ulna, the true value of the bending stiffness can be 

deterained. I coaparison of this value with the foraer values 

will reveal whether or not the aisalignaent is the only cause of 

the rotational resistance at the support, and which positioning 

will ainiaize or elininate it. Several sets of such experiaents 

will aid in establishing a standard, systeaatic aethod of 

* 

positioning for all future in vivo aonkey foreara tests. 

The in vitro experiaents suggested in this section, as well 
as those discussed in Section VI. 1 should also be perforaed on 
aonkey tibiae, huaan cadaver ulnae and any other type of 
speciaen to be routinely tested £n vivo , llthongh the aodeling 
concepts applied to the foreara of a aonkey are also applicable 
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to monkey legs and human foreans, the geoaetry of the supports 
in each case is guite different. A standard, systenatic aethod 
of positioning is also needed in these cases. 

3. Suggested Bonifications to the Test Procednre 

The iapedance aeasuring procedure currently being used at 
laes Research Center has one aajor flaw: the support conditions 
of the specinen are not being controlled. Since the DPBI is just 
as sensitive to the support conditions as it is to the bending 
stiffness of the bone, the support conditions aust be known in 
order to deternine the bending stiffness. If the boundary 
conditions of the aatheaatical aodel are not a good 
representation of the support conditions of the physical system, 
then the value obtained for the bending stiffness will be in 
error, possibly as aocb as an order of magnitude. 

Two aodif ications to the inpedance aeasuring procedure are 
reconoendei . First, the positioning procedure to be established 
by the experiments suggested above should be adopted as part of 
the procedure for each DPBI test. This will reduce, if not 
completely eliminate the resisitance to rotation at the 
supports. Second, Thompson's procedure, involving the use of 
plaster pads under the wrist and elbov should be readopted. This 
will eliminate the translation allowed by the putty at the 
supports ' (Section 711. B) . The result of adopting these 
aodif ications is that the support conditions will be 
sufficiently controlled to obtain repeatable accurate results. 

One further recomaendation which aay prevent the production 
of aeaningless DPBI data is suggested. A standard, such as the 
aluminum beaa (Section 71. B) , should be used to insure that the 


impedance measuring equipment is operating properly over the 
frequency range of the test. Each tine DPBI tests are conducted, 
the DPBI of the standard should be Measured and the data briefly 
examined. For example, using the aluminum bean shovn in Figure 
6.19, vith the pins lubricated vith a light oil, the general 
shape of the DPB1 data should be as shovn in Figure 6.20. The 
■iniaua point should occur at approxiaately 450 Hz and the 
aaxiauB point at approxiaately 2800 Bz. The static stiffness 
should be 5.35 x 10 T dyne/ca which corresponds to a DPSI of 8.5 
x 10* dyne s/ca at 100 Hz. If these specifications are not set 
to within a few percent, then the iapedance aeasuring eguipaent 
should be further checked for aalf unctioning. 

4. Concluding Remarks 

The impedance aeasuring procedure developed by Thompson 
(Section l.E), vith recoaaended modifications discussed above, 
can be used to generate an accurate, repeatable set of DPHI data 
for a foreara or leg. A systeaatic, user oriented analysis 
procedure has been developed and prograamed on a digital 
computer. The coaputer program, listed in Appendix F, employs 
the aatheaatical model, developed in Chapters III and IV, and 
the SID&, developed in Chapter V. The aatheaatical model 
consists of a unifora, linear, visco-elastic, simply-supported 
Euler-Bernoulli beam to represent the bone; a visco-elastic 
foundation with mass to represent the tissue surrounding the 
bone; and a spring between the bean and driving force to 
represent the skin between the bone and probe. The SIDA 
determines values for the aatheaatical model which best simulate 
the DPBI data using an iteration scheme to minimize an error 
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function. The error function is sinilar to that which is used in 
a classical least-squares curve fit. Due to the resemblance 
between the mathematical model and the physical system, the 
parametric values which produce a good simulation of the DPHI 
will infer the material and geometrical properties of the 
physical system. 

One of these properties, the bending stiffness of the bone, 
was shown to correlate quite well with its breaking strength, at 
least for normal bones (Borders, Petersen and Orne, 1977; Jurist 
and Poltz, 1977). Breaking strength is a good measure of bone 
integrity and therefore may be a good indicator for many bone 
disorders such as osteoporosis. However, more correlation 
studies are needed to determine the effects of various bone 
disorders on the stiffness and strength of bones. 

Bone mineral content (BHC) is currently being used in 
ongoing experiments to monitor changes in the bones of monkeys 
during prolonged hypodynanic restraint (Young and Tremor, 1976). 
Impedance testing is the only feasible technique currently 
available as a possible countermeasure to BHC. The impedance 
measuring and analysis procedures presented here can be used in 
conjunction with measurements of BBC to better define the 
condition of the bone being examined. 

Young and Tremor (1978) report an average of 3.5 percent 

* 

loss in femoral BHC in ten restrained monkeys over the 
relatively short time period of one month. Vbedon et al. (1976) 
reports changes in BHC of 7.9 percent in the os calcis of 
astronauts after 84 days in a weightless enviroment, in spite of 
a rigorous exercise program. These changes are significant 
although they occurred during a relatively short period of time. 


Such larger changes are expected to occur over longer periods of 
weightlessness, e.g., daring a 1.5 to 3 year trip to Bars, or in 
a severe case of bone disease such as osteoporosis. 

Although the percent changes in bending stiffness vhich 
occur with various bone disorders have not been aeasured, they 
are expected to be at least as great as those found in BBC. 
Bending stiffness is proportional to the fourth order of the 
cross sectional dinensions while BBC is proportional only to the 
second order, i.e., 

El * B C, d« BBC = BBD A « BAD C*dz (8.1) 

where BAD is the bone nineral density, 

A is the area of the cross section, 
d is a cross sectional diaension, 
and c, , c 2 are constants of proportionality. 

Therefore, the bending stiffness is aore sensitive than the BBC 
is to changes in geoaetry. If percent changes in aodulus of 
elasticity are of the sane order of aagnitude as percent changes 
in BBD, then bending stiffness will actually be a aore sensitive 
indicator than BBC. Bence, the expected percent changes in 
bending stiffness are greater than those cited above for BBC and 
greater yet for aore severe cases. Bith the recoaaendations 
discussed above taken into account, the iapedance aeasuring 
procedure is accurate and repeatable enough to detect and 
aeasure these changes. 

A technician in the clinical setting, can carry oat the 
iapedance testing procedure and run the coaputer prograa to 
deteraine the bending stiffness of a bone and interpret the 
result in teras of a particular bone disorder, all with a 
ainiaua of training. The test takes only a few ainutes and is 
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entirely noninvasive . Two developments are needed to ascertain 
the feasibility of this technique. They are: (1) to develop a 
systeaatic positioning procedure, and (2) to develop the 
correlations between BflC, bending stiffness and various bone 
disorders. Both of these are quite achievable. 


CHAPTER IX 


APPENDIX 


JU I8PEDANCE eopations 

The driving-point nechanical impedance (DPHI) of a single- 
degree-of-f reedoa oscillator is 
Z* = c ♦ i (ap - K/p) 

The DPMI of a bean is of the forn 
Z* = 2E*IA*/[ ipf (AL) ] 
where 

A* = /i*p 2 /E*I 

and f (AL) is a function which depends on the boundary 
conditions. For each set of boundary conditions listed in Table 
3.1, f (AL) is as follows: 

1. Sinply-supported 

f (AL) * sinAa sinAb/sinAL - sinhAa sinh Ab/sinh AL 

2. Rotational spring on one end 

f (AL) * [ (sinAa ♦ x) (sinAb sinhAL * (!) 

f 

- (sinhAa ♦ a) (sinhAb sinAL ♦ p) ]/ 

i 

(sinAL sinhAL ♦ Y) 

where * * k, (coshAa - cosAa)/2E*I^ 

(? = k, (sinAb coshAL - sinhAb cosAL)/2E*IA 
y = k, (sinAL coshAL - sinhAL cosAL)/2E*IA 




1 



o 
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3. Rotational spring on each end 

f (AL) = [ (sinAa ♦ a) (sinAb ♦ (3) (sinhAL ♦ Y ♦ £,) 

-(sinhAa ♦ a) (sinhAb ♦ 0) (sinAL - V - £ t ) 

-(sinAb ♦ p) (sinhAa ♦ ft) (V ♦ S 4 ) 

-(sinAa ♦ a) (sinhAb ♦ f) (y ♦ &,) ]/ 

[ (sinAL - V - £, ) (sinhAL ♦ Y ♦ s t ) ♦ (V ♦ S 5 ) (/ ♦ SJ ] 
where a 8 k, (coshAa - cosAa)/2E*IA 
0 8 k 2 (coshAb - cosAb)/2E*IA 

Y 8 k, k x (sinhAL ♦ sin AL) / (2E*IA) 2 

S, * (k, ♦ k 2 ) cosAL/2E*IA 

Si 8 (k, ♦ k 2 ) coshAL/2E*I A 

S s « (k, cosAL ♦ k^ coshAL) /2E*IA 

£ 8 (k* cosAL ♦ k, coshAL) /2E*I A 

Q. Translational spring on one end 

f (AL) = sinAa sinAb/sinAL - sinhAa sinh Ab/sinhAL - p 2 /S 
where $ = sinAb/sinAL ♦ sinhAb/sinhAL 

S = coshAL/sinhAL - cosAL/sinAL - 2k,/E*IA* 

5. Translational spring on each end 

f(AL) * sinAa sinAb/sinAL - sinhAa sinh Ab/sinbAL 
- * K*6, ♦ 2 apy)/U,S t -V*) 

where m = sinAa/sinAL ♦ sinh Aa/sinhAL 
f? 8 sinAb/sinAL ♦ sinhAb/sinhAL 

Y 8 1/sinhAL - 1/sin AL 

S, 8 coshAL/sinhAL - cosAL/sinAL - 2k,/E*IA* 

8 coshAL/sinhAL - cosAL/sinAL - 2k z /E*IA 3 
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6. Translational spring on an extended bean 


f (AL) * sinAa sinAb/sinAL - sinbAa sinhAb/sinbAL 

-p/[ IV ♦ *£)/(£ - *) - SJ 

vhere p * sinbAb/sinbAL - sinAb/sinAL 

Y * 2(cosAe cosbAe ♦ 1)/sinAe sinble 
* * 2k 3 /E*IA* 

c * cosbAe/sinble - cosAe/sinAe 
6 * cosbAL/sinhAL - cosAL/sinAL 
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The transfer equation across the driving force is 
{!<*♦>) = {Y< r >} ♦ {F} 



and the (♦) and (-) superscripts refer to the state variables 
just to the right and left of the driving force applied at the 
kth node. Let 

* 

[0] = [ T* ] [T k „ ]...[!*] [T, ] 

P] = [T„] [T„ H T k , l ] [T k *,] 

C S3 * [V] [0] 

then the following two satriz equations are obtained by 
successive substitutions froa one transfer aatriz egnation to 
the nezt 

(Y< k ->) * [0] (I e J 

= [S] Ho) ♦ [V] (F) 

These tvo aatriz equations represent eight algebraic equations 
of twelve state variables. Four of the state variables sust be 
known froa the boundary conditions leaving eight unknown state 
variables. 

Any set of classical or non-classical boundary conditions 
can be applied to these eight equations. The simply-supported 
boundary condition states that 

y, - 7 n - 0 b, - a„ « 0 • 

After applying these, the first, fifth and seventh equations are 

Jk c °'1 0 o ♦ 0 u V o 
0 = s a e c ♦ s„v e ♦ ?„F 
0 B S^0 a v S u v 0 ♦ IjjF 

The solution for y K , after eliainating 0 O and T 0 froa these 
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three equations is 

y K ■ F [O a (V, 4 S, 4 - ? U S^) ♦ 0 M (V )4 S,t- V^s^) l/(S (i S 34 - S (4 S,J 
Finally, the DPHI is 

2* * F/ipy* 
or 

2* * (“i/p) ~ S l4 SjJ /[0|i(1}»S l4 “ 7|*S 34 ) ♦ D |4 (V /4 Sjj. “ ^ 4 S| t ) ] 

Since the exact solution for each elesent was used, the 
accuracy of the total solution is as good as the accuracy of the 
step function approxination of the taper. 
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£t DEPENDENCE OP THE FOONDATION PARABETEBS ON THE BINIHPB 
POINT OF AN I H PEDA NCE PLOT 

The ainiaua driving-point aechanical iapedance (DPBI) is to 
be deterained for various coabinatiODs of the values of the Bass 
per unit length p { , and the daaping ratio X;$ of the foundation. 
Several DPBI plots, siailar to those of Figures 4.13 through 
4.16 were generated. 

DPBI plots are generated by evaluating the DPBI equations 
at a finite nuaber of points and joining those points with a 
sequence of straight line segaents. 1 large enough nuaber of 
points are taken to give the DPBI plots the appearance of snooth 
curves. The values of the DPBI and the forcing frequency for 
each point are listed in the computer printout associated with 

each plot. The true ainiaun point aay not occur precisely at one 

of these points. In such a case, the true ainiaua point occurs 
at sose frequency betveen the frequencies of the lowest DPBI 
listed and an adjacent point. The true ainiaua DPBI is lover 
than either of these points. See Figure C. 1. A good 

approxination to the true ainiaua DPBI is obtained froa the 

values of the DPBI and forcing frequency of the lowest point 

listed and its two adjacent points as follows: 

/ 

(x 0 ,y 0 ) describe the coordinates of the true ainiaua 

f 

point, i.e. , 

Xo * Jo * 

Siailarly , let (z £ »y 2 )» (*,#!, ) end <x,,y 3 ) describe the 

coordinates of the lowest listed point and its two adjacent 
points, respectively. See Figure C. 1. Approximate the DPBI 
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equation by a quadratic equation in the region around these 
points 

y ■ A ♦ Bx ♦ Cx* 

1, B and C are constants which can easily be found by solving 

y, ■ A ♦ Bx, ♦ Cx* 

y* * A ♦ Bx* ♦ Cx| 

jj « i ♦ Bx, ♦ CxJ 

The ainiaua point frequency is found by seting the derivative of 
the DPMI equations equal to zero 
y* * B ♦ 2Cx * 0 
X 0 * -B/2C 

The ainiaua DPMI is obtained by replacing x with the expression 
for x 0 

y D * 1 - B7«C 

The ainiaua DPMI vas deterained in this way for several 
DPMI plots, each generated with a different coabination of 
values of p f and fy. The results are tabulated in Table C.l. 

Por the case where p* * 0 and T/ * ® (i.e., no foundation), 
the ainiaua DPMI is given by 
Z M(U u,/K - 2 JT 

(see Sections IV. A and 17. B) . Therefore a reasonable fora to 
assuae for the ainiaua DPMI is 

ft 

Z„„u>/K « 2r ♦ f (p M Y f ) 

where f (p f * ? f ) is a function of pj and whose value is zero at 
P; * 0 and = 0. Values of this function are found froa values 
of the ainiaua DPMI by subtracting 2£. The results are tabulated 
in Table C.2. 

Bote that the values of f(p*#? ; ) seen to increase linearly 
with P//p* Assuae that the dependence of p; on f (p; , fy) is in 
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fact ’ii.car, i.e., assuae 

Z miu w/K « 2f ♦ p,/p g (r#> 

where g(^) is a fonction of f f whose walae is zero at ^ *0. 
Values of this function are found froa walues of f(p (t?f) 
dividing by their respective values of p f / p. The results are 
tabulated in Table C.3. 

To deteraine the fora of the function g(^), its values 
were plotted on a log-log grid, ill points were found to lie 
very close to a straight line. Therefore g(^) is of the fora of 
a power of f t§ i.e., 

. Z m *w/K * 2? ♦ p,/p A t? 

where n is one-third for a fixed foundation and one-half for a 
free foundation. To find the values of A, a least-sguares 
technigue was c-aployed. The best values for A were found to be 
one-fourth for a fixed foundation and three-fourths for a free 
foundation. 

The relationships which approxiaate the dependence of p f 
and on the ainiaua DPfll are 
Zm,v^/K * 2? ♦ V« Pf/P 

and 

Z^w/K « 2J ♦ 3/a p,/p 

for a fixed and free foundation, respectively. These 
relationships can be used to deteraine approxiaate values for 
the tissue- parameters of a vibrating foreara or leg systea 
directly froa its DPfll data plot. Such an approxiaation is 
necessary to establish the initial guess for the systeas 
identification algorithn discussed in Section V.C. 


Di THE BIHIHOM AMD HAIIHOH POINTS OF AH IMPEDANCE PLOT 


Expressions for the ainiaua end aaxiaua points of the 
driving-point aechanical iapedance (DPMI) plot are to be found. 
The lengthy analysis will be outlined briefly here. 

The DP HI of a single-degree-of-freedoa oscillator in series 
uith a spring (see Figure 4.21) with % # 0 is 
2* = [1/(8ip ♦ c ♦ K/ip) ♦ 1/(k/ip) J~ l 
After replacing • by K/w 2 , c by 2KJf/o and perforaiug several 
steps of algebra, this equation becoaes 

2*S(0 2 -1) -2£S0(02-i-s) 

-i[«?*sp*S(p*-1) (02-1-S) ] 

2* = K/o, 

(f»2-1-S) 2+CJ20 


where S = k/K and 0 « p/uj. The aagnitude of the DPHI is 


[2£S(02-1) -2£S0(0=-1-S) ]2 
\ <-[4j s 2SfS+S(02-1) (02-1-S) ]2 

2 » K/ ~ 

0(£*-1-S)Z + 42f20 

To find the ainiaua and aaxiaua points, take the derivative 

of the aagnitude of the DPHI and set it equal to zero. 

((3(02-1-5) 2 *4r20) 

X {[ 4£S0-2£S (02-1-S) -4^50*3 
x[2?S(02-1)-2rS£(/S*-1-S) ] 

♦[ 4^S + 2S0(2|52-2-S) 3 
X[4^^S0*S (02-1) (02-1-S) 3) 

-{(^2-1-S) 2 + 402 (02-1-S) + 4/2) 

* l[2^S (0 2 -1) - 2 0(02-i-s) 3 2 
+ [4^2Sp+S(02-1) (02-1-S) 3 2 ) 
d2/dp « K/ * * 0 


(0 f£2-i- S ) *±am 

f 2/S(p2-i)-2^S0(£2-1-S) 3 2 
xli+th^s^+s^z-ii (0 2— i — s) 32 

The denoainator is positive and therefore ton-zero for all 
positive values of {$, y and S. Therefore the nuaerator aust be 
set equal to zero. The expression in the nuaerator, when 
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■ultiplied oat, is a sixth order polynonial in ?*. As an 
alternative to the difficalt task of solving it, Taylor series 
expansions of p 2 with respect to can be found which satisfy 
the sixth order polynonial equation. 

Assune that the solutions for p 2 exist and are of the forn 

■ S ♦ 1 » £ a.?*" 

PJ • 1 ♦ t b "?'" 

n-i 

where = P*a,/co and * P M . W /^- These equations will produce 
the correct solutions for * 0 according to equations (4.23) 
and (4.27). Substitute the assuned fora of the solutions into 
the nunerator of the equation. The coefficients of the constant 
tern and the and terns are each set equal to zero. In each 
case, the constant tern was found to be identically equal to 
zero, indicating that equations (4.23) and (4.27) are actually 
the correct first order approxinations to the solutions. The 
equations obtained fro* the and *£* terns are solved to obtain 
the first two unknown coefficients of each of the Taylor series. 
Hence, the first three terns of each of the Taylor series are 
found to be 

= S ♦ 1 ♦ 2/S (2+S)/(HS) t* 

- 2/S* (2 + S)/(US®) (4v16S*13S 2 + 4S*) •£* ♦ ••• 
ft** - 1 - VS t* ♦ 8/S® (2+3S) - ... 
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Ei DERIVATIVES OF THE I HPBDANCE 

The driving-point mechanical impedance (DPHI) equation to 
be differentiated is 

Z* = {ipf (AL) /[ 2EI (1 ♦ 2i~$p/uj) A 3 ] ♦ ipA} -1 
where 

f (AL) = sinAa sinAb/sinAL - sinhAa sinhAb/sinhAL 
% * [ttVL 4 p 2 /w* ♦ p2p^/EI g (V^) ]*A (1 ♦ 2i/p/u^-«A 
r - pn/aJf / -1 1 «■ 2iJJp 

and the function g (f ) , depends on the type of foundation 
included in the model. Three cases are considered. Case A: no 
foundation. The function g (¥>) , is zero and A recuces to 
A = t r/L (p/uj) ^ (i ♦ 2i?p/u>)-iM 
Case B: fixed foundation 
gCV) = —cot V / y 
Case C: free foundation 
g(V) = tanr/2 / f/2 

Define X and T as the real and iaaginary parts of the 
inverse of the conplex DPHI, respectively, i.e., 

Z* * (X ♦ il) **» 

The magnitude of the DPHI is 
Z * (X 2 ♦ T*)-V2 

r 

The derivative of the lagnitude of the DPHI with respect to one 
of the nodel parameters is 

dZ/dP «=-(!*♦ I*)-** (X dX/dP ♦ I dl/dP) 
or dZ/dP * -Z» (X dl/dP ♦ I dl/dP) 

where P represents any one of the nodel parameters. The value of 








X, T and their derivatives are calculated froa 13 

X = Peal (1/Z*) dl/dP = Beal[ d (1/Z*) /dP] 

I = laag (1/Z*) dl/dP « Iaag[d (1/Z*)/dP] 

Since the DPHI is a function of EI r to, k and A; and A is 
a function of El, co, Jf, p, , 10 ^ and £ f ; the derivatives are of 
the fora 

d (1/Z*) /dEI = 0(1/Z*)/dA dA/dEI ♦ d(1/Z*)/dEI 
d (1/Z*) /dio = 9(1/Z*)/dA dA/du> ♦ d(1/Z*)/dcu 
d (1/Z*) /d£ = d (1/Z*) /dA dA/djT ♦ OO/Z*)/^ 
d (1/Z*) /dp| * d(1/Z*)/0A dA/d|i* 
d(1/Z*)/du; f - y ; 1/Z*) /dA dA/duJ f 
d (1/Z*) /d^ = d (1/Z*) /dA dA/d£ f 
d (1/Z*) /dk ■ a (1/Z*) /Ok 
The partial derivatives are 

O(1/Z*)/OEI = -ipf(AL) [ 2 (El) * (1 ♦ 2iyp/w)A 3 ]-* 
d (1/Z*) /Ou) = -p2f(AL) £uj 7 El ( 1 ♦ 2i*p/u>) 2 A 3 ]-» 
d (1/Z*) /3£ = p 2 f(*L) [wEI(1 ♦ 2i?p/u;) 2 A 3 3-» 

0 (1/Z*) /3k * -ip/k 2 

d(1/Z*)/3A *= -3ipf(AL) £ 2EX (1 ♦ 2i£p/w) A* ]"> 

♦ ip df/dA [ 2EI (1 ♦ 2i£p/u>) A 3 1 

where 


13 Since X and T are real continuous functions, and i = - 1 is a 

constant, ' the distributive property of the derivative holds, 

1* 6a f 

d (X ♦ iT)/dP * dl/dP ♦ i dl/dP 
Bence Beal[d(X ♦ iI)/dP] * dl/dP 
and Iaag[d(X ♦ iT)/dP] «= dl/dP 


1 




o 
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Df/dA = [a cosAa sinAb ♦ b sinAa cosAb 

- L sinAa sinAb cosAL/sinAL] / sinAL 
- [a coshAa sinhAb ♦ b sinhAa cosbAb 

- L sinhAa sinhAb coshAL/sinh Al] / sinhAL 

Since /l is a function of El, u/, p^ and Y; and V is a 
function of u)± and the derivatives of A sure of the fora 
dA/dEI * 3 A/3 El dA/dpj « 3 A/Op, 

dA/dw = 3A/3u> dA/dw f * 3A/3V dr/d to, 

dA/d£ = 3A/3? dA/d& = dX/df d¥yd* f 

The partial derivatives of are 
Case & 

3A/3EI * 0 

©A/Ou) = -rr/2uiL (p/o) (1 ♦ 2i^p/u>) “ l '* 4 

♦ iri*p/2Lu* (p /u>)** (1 ♦ 2i]?pAo)-** 

©A/©£ * - jrip/2wL (p/u>) fc* {1 ♦ 2i£p/w)-** 

Cases B and C 

OA/OEI = -p*p f /4(EI)* g(^) (1 ♦ 2i^p/u>) “tA 

[rr 4 /L 4 p 2 /oo* ♦ p 2 Pf /El g(t) }"** 

©A/0^ = -1/2 to rr^/L 4 p 2 /co 2 (1 ♦ 2i2pA>)”** 

[rr 4 /L 4 p 2 /u* ♦ P 2 p f /El gMJ-** 

♦ i£p/2to* (1 ♦ 2i£p/io)“** 

£t t */ i * p 2 /w* ♦ p 2 » f /Bi g(t) ]•* 

3A/35 * -ip/2 to (1 ♦ 2i$p/to)-** 

[rr 4 /L 4 p 2 /a* ♦ p 2 p*/EI gW) ]»* 

3A/3p f « p 2 /4EI g(f) (1 ♦ 2i*p/io)~ l/4 

[tt 4 /i 4 p 2 /u»* ♦ p 2 p,/ei g{f) r w 

©A/Of = p*p,/4EI dg/d? (1 ♦ 2i£p/aj)-t* 

[tt 4 /l 4 p*a* ♦ p ? p*/ei g(r)] - ** 

where 




dg/dy = [V esc 2 '*' ♦ cotf]/V 2 
dg/d'f' = 1/2 [V/2 sec*Y72 - taey/2 ]/ (Y/2) 2 
The derivatives of y are 
dv/du;* * -1/ujf Ptt/o)| (1 ♦ 2i^p/u>4 

-i^p/u^ 2 pTT/^ (1 ♦ 2i$ f p/a^)-V* 
dV'/dff * -ip/uy* pTr/^ (1 ♦ 2if f p/o^)-** 


(case B) 
(case C) 
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Zs. THE COHPOTER PROGRAM 

I listing of the Fortran computer program which determines 
the parametric values of the mathematical model used to simulate 
a set of driving-point mechanical impedance data from a forearm 
or leg vibration test is given. Ill of the function subroutines 
required by the program are not available in double precision. 
Therefore, five function subroutines have been written to 
accommodate the main program. They are also listed. The 
subroutine DGELG from the IBS Scientific Subroutine Package 
(*SSP) is used to solve the system of linear algebraic equations 
within each iteration of the systems identification algorithm. & 
listing of DGELG can be found in IBfl (1968). 
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JHIS_PEOGEAH .XBI10IS .AN.-IT.ER.ATIl£_PitQCEDOBE JT OJT ON Vj;.P5j:_01L_IilE 
CORRECT VALUES OF SHE PARAMETERS IN A VIBRATING LONG DONE ' ' 


EXPERIMENT, BY MINIMIZING THE PERCENTAGE ERROfi IN THE MAGNITUDE 
OF THE IMPEDANCE. 


THE INPUT DATA BUST BE ARRANGED AS FOLLOWS: 


15A4 


CA RD 1 T I TLE 15A» 

CARD 2 IE NGTH AND LENGTH- TO -PRO BE- 

LOCATION RATIO FREE 

CARD 3 BOON DARI CONDITION OF TISSUE 12 

CAED 4 NUMBER OF DATA CARDS TO FOLLOW 12 

THE REST OF THE CARDS CONTAIN THE FREQUENCY AND THE HAGNITODE AND 
PHASE ANGLE OF THE IMPEDANCE, ONE POINT PER CARD, IN FREE FORMAT. 


CAE D 

CARD 


FREE 

12 

12 


THE SIX PARAMETERS IN THIS MODEL ARE: 


BEI 


BUN 

TMU 

TUN 

I ZETA 

K 

BZET fc, TH E 
VALUE. 


ST I F FNESS OF THE BONE 

NATURAL FREQUENCY OF THE BONE 
HASS PER UNIT LENGTH OF THE TISSUE 

NATURAL FREQUENCY OF THE TISSUE 

LAMPING RATIO OF THE TISSUE 
STIFPNESS OF THE SKIN 

DAMPING RATIO OF THE BONE, IS HELD AT A CONSTANT 


THE -FOUNDATION IN THE MODEL, WHICH REPRESENTS THE TISSUE, CAN 
HAVE EITHER A FIXED OR FREE BOUNDARY DEPENDING ON THE VALUE ON- 
CARD 3. -1 CORRESPONDS TO A FIXED BOUNDARY. 1 CORRESPONDS TO A 

FREE BOUNDARY. 


THIS PROGRAM CONTAINS ROUTINES WHICH 'LOOK* AT THE DATA AMD 
CHOOSE INITIAL SETS OF PARAMETER VALUES. 


THE ITERATIONS ARE CARRIED OUT IN THREE PHASES: 

1. A FOUR PARAMETER MODEL IS EMPLOYED TO OBTAIN A GOOD 
’APPROXIMATION TC THE BONE AND SKIN PARAMETERS. 

2. THESE ARE HELD FIXED WHILE A GOOD APPROXIMATION TO THE 
TISSUE PARAMETERS IS OBTAINED FOB A SIX PARAMETER MODEL. 

3. ALL SIX PARAMETERS ARE ALLOWED TO VARY TO OBTAIN THE FINAL 
SET OF PARAMETERS FOB THE SIX PARAMETER MODEL. 


DECLARATION STATEMENTS. 


COMPLEX* 16 DCMPLX,CDSQBT,CDTAN,CDS INU,CDSIN ,CDABS,CDCOSH ,CDCOS 
COMPLE X*~16 DZI (6) 

COMPLE I* 1 6 AEG ,LAM PA ,ZA ,ZB ,ZL, ZBI, ZTI , ZC,BQ,TQ, Q, 

ZBI 1, Z BI2, DZI Di, COT, CSCS, SECS 
R EAL*B D B L 2 , D E E A L , DIM AG, DATA'S, DABS 

REAL*8 8 (60 ) , P (60) ,22(6 0) ,PHIE(60) ,Z(60) ,PUI(60) ,DZ (60,6) ,DP(6), 

A 6 (6,6) ,A4 {4, 4) ,A3 (3,3) ,B (6) ,DP6 (6) , DP4 (4) ,DP3 (3) ,DX (6) ,DY (6) 

BEAL*S PI, Bl, B FAT 10, BA, Bfa, BEI, BUN, BPN ,'LZETA, 

TflU,TUN,T?N ,TZETA,K,KK,ZMIN,ZMAI,WHIN,WMAX,X, Y, ERROR, ERROLD, 
SBEI,SfiUN,SBZETA,SK,STMU,STUN,STZETA,& 

B E A L * 3 F M I fl J / 1 * 1 / 

INTEGER TITLE(IS) < rlGfMdr 

BEAL W P (60 ) ,ZP (60) , PHIP (60) , ZEP (60) ,PHIEP(60) qp ^AGfi IS 
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x)- 


BEAD IN DATA. 
PBINT74 


PI S 3. 1 4 1 59D 0 
BEAD (5, 5) TITLE 

PO E 1ATJ15A4J 

BEAD (5,FHT) EL , BE ATI 0 
BL=DABS (BL) 

J ZJL B SATIO.LE. O.DO.OE.B5ATIO. GE. 1.D0) GO TO 6 
GO TO 8 
PBINT7 

POEM AT I* THE TALOE GIVEN TO THE LENGTH-1 


THE TAL DE GIVEN TO THE LENGTH-TO-PB ODE-LOCATION*, 
BATIO RUST BE BETWEEN ZEBO' AND ONE.*/) 


STOP 

CONTIROE 


BA-BL* BBATIO 
BB=BL-BA 

BEAD (5, 20) IBC 

IF (I3C. E Q. 0 ) GO TO 16 
IBC=ISJ.GN (1 ,IBC) 

I P ( I B C^E^. -1) P R INTI 2 


FOBHAT (' 


IP (IBC. EQ. 1JPRINT13 


THE BOUNDARY OF THE POONDATION IS FIXED (ULNA) • */) 


THE BOUNDARY OF THE FOUNDATION IS FBEE (TIBIA) 


THE FOUNDATION IS NOT INCLUDED IN THE HODEL* 


FO BflATl! THE BOUNDARY OF 1 

GO TO 18 

PBINT17 

__PO R N A Tj^ TH E FO'JN DAT ION IS 

CONTINUE ~ 

££AD(5,20) N 

PO B SAT (I 2J 

" IF (N.LT. 8) GO TO 21 
IF (N . GT. 60) GO TO 23 
GO TO 25 
PRINT22 

FOBHAT (' A RINIBUH OP EIGi 

STOP 

PB IN 12 4 

FOBHAT (» A HAIIMUH OP Sill 

STOP 

CONTINUE 
DO 26 1=1,1 

B£AD(5,PHT) B (I) , Z E ( I) , PHI E (I) 
H (I) = DABS (B (I) ) 

ZE(I) =DABS (ZE (I)) 

P(I) =W (I)*2.D0*PI 
PBINI7 4 


A RINIHUH OP EIGHT DATA POINTS IS BEQUIBED. •/) 


A HAIIMUH OP SIITI DATA POIHTS IS BEQUIBED.*/) 


PHASE 1 


DETEBHINE INITIAL SET OP PAB1HETEBS. 

Tk=o.dq 

DO 27 1=1,4 
KK=KK*ZE (I) *P(I) 
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_FK = KK/4.DQ _ 

sW 


NN=N/2 

>8 


WMAX=R (NN-1 ) 

>9 


ZMAX=Z E (NN- 1) 

>0 


DO 30 I=NN,S 

il 


IF (Z E ( I) . GT.ZB AX) GO TO 29 

>2 


GO TO 30 

>3 

29 

ZMAX=ZE (I) 

>4 


BMAI = B (I) 

>5 

30 

CONTINUE 

> 6 


HN=N/4 

>7 


HNN=3*N/4 

.8 


WMIN=W (NN-1) 

.9 


ZMIN=ZE(NN-1) 

'0 


DO 32 1= NN, NNN 

’1 


IF (ZE (I). LT.ZHIN)GO TO 31 

'2 


GO TO 32 

3 

31 

ZNIN=ZE(I) 

'4 


WMIN=H (I) 

5 

32 

CONTINUE 

6 


K=KK* (BMAX/BHIN)**2 

7 


KK=1.D0/( (1.D0/KK) -(1.D0/K) ) 


8 BEI= (BA*BB) **2/3. D0/BL*KK 

9 BW N= NHIN 

0 B Z ET A =2H IN»PI»BMM/K K 

1 P2IKT33 ‘ 

2 (j 33 FORMAT (/• THE INITIAL SET OF PARAMETEBS IS:*/) 

- PRINT77 

PR INT7 8, BLI ,B V N ,BZE1 A, K 

THIS IS THE BEGINNING OF THE OUTSIDE LOOP, EACH BUN THROUGHTfilS 
LOOP CONSTITUTES ONE ITERATION, 

EBRO B= 1. D20 

6 HM=0 ' 

7 34 HM=MH* 1 


4 

C 

C 

c 

c 

5 


3 

9 

0 

1 

2 

3 

» 


j 

7 


3 r' 


c 

c 

c 
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c 

c 


CHECK EACH PARAMETER FOB THE NON-NEGATIVITY CONDITION. 

I P (BEI. LT. 0. D0) GO T O 35 

IF ( B k N . LT. 0. DO) GO TO 35 
IF (BZETA.LT. O.DO)GO TO 35 

IF(K.LT. 0.D0)GQ TO 35 

GO TO 37 
HMM= HH— 1 

PRIWT36, HHM ‘ 

FORMAT 1/51 , • A NEGATIVE VALUE BAS OBTAINED FOR ONE OR MORE * /5I, 
•OF THE PARAMETERS ON ITERATION NUMBER • # I2, , .*/5X, 

•THE CURRENT PARAMETER VALDES ABE:*/) 

PRINT?' 7 

paiNT78,BEI,BSN # BZETA,I 

6°. T0 -_ 4 .® i^lOXAL PAGE fci 

CONTINUE OF POOR QUALITY 

CALCULATE Z AT EACB FREQUENCY. 
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9 

0 

1 

2 

3 

4 

5 


7 



3 

) 

) 

1 


J 

\ 

j 

> 

1 

r 


i 

i 

» 

1 


~o 


c « 

BPN=BWN*2.D0*PI 
DO 45 1=1, H 

BQ=DCH ?LX (1.D0 , 2 . 00 »B 2 ETA* P ( I) /B PH) ' 

LARDA= DCRPLX (PI/BL,0. DO) *CDSQET (DCRPLX (P (I)/BPN,O.DO) /CDSQET (BQ) ) 
ZA= LAM DA * DCS PL X (BA,O.DO) 

Z B = I. A H D A_*D CMP L X (BB ,0 « DO) 

ZL=LAH DA *DC KPLX (EL,Q. DO) 

ZBI1=DCS PLX (0. DO,-.5DO*P(I) /BEI) /LAH DA**3/BQ 
ZB 1 2 s CDS IN R (Z A) *CD SI NR (2B) /CDSINHfZL) 

1 -CDSIN (2 A) * CD w IN (ZB) /CDSIH (ZL) 

ZBI=ZBI1 *ZBI2 

ZTI=DCRPLX (0. DO # P(I) /K) 

ZC=ZTi*ZBI 
X=DBEAL (ZC) 

If D IRA G (ZC) 

ZC=DCRELX(1.D0,0. DO) /ZC 
Z (I) s CDABS (ZC) 

C 

C CALCULATE THE CEBIV1TITES OF Z AT EACH FBEQUENCI* 

C 

DZIDL= ( DC H P LX ( E A 1 0 ,_D0 )_*C DCOSH (2 A) *CD SIH H_( Z B] 

2 ♦DCRPLX (DU, 6 .DO) • CDSXNIi'(ZA) *CDCOSH (ZB) 

3 -DCRPLX (BL,0.D0) ♦CDSIKh (ZA)*CDSINH (ZB) 

4 ♦CDCOSH (ZLJ/£DSINH (ZL) ) /C DSIK H (ZL) 

DZIDL= DZIDL- (DCRPL X ( EA, 0. D0)'*CbCOS (ZA) *CDSIN (ZB) 

2 ♦DCMPLX (BB,0- DO) *CDSIN (ZA) *CDCOS (ZB) 

3 -DCBPL X ( 3_L, 0. DO) ♦ CDSIN ( Z_A)_* C DS IHJ ZB) 

4 *CDCUS (ZL) /CDSIK (ZL) ) /CDSIH (ZL) 

DZIDL-DZ IDL*Z B1 1+ DCMPLX (-3. DO, O.DO) /LARD A*ZBI 

Q= DZ I DL* DC M PL X (-.2 5D 0 ,0.D3j ♦LAMDA/BQ-ZBI/BQ 

DZ 1(1) =- ZB I /DC 2 PL X (B EI , 0 . D 0 ) 

DZI (2) =Q* DCRPLX (0. DO, -2. DO*BZETA*P (I)/BPH**2) 

1 ♦DZIDL*DCMPLX (-. 5D0/EPH , 0. DO) *L ANDA 
DZI (3) =Q*bC fl?L X (0. DO , 2 . DO * P (I) /B PH ) 

DZI (4) =-ZTI /DCRPLX (K,O.DO) 

DZI ( 2) = DZI (2) ♦ DCRPLX (2. DO*PI,O.DO) 

DO '45 J= 1 , 4 
DX (J)=DBEAL (DZI (J) ) 

Dl (J) =DIRAG (DZI (J) \ 

45 DZ (I, J) =-2 (l) *»3* (EX (J) *X ♦ DY (J) * X) 

C 

C CALCDL ATE A HD PEIHT THE EBBOB POHCTICH, 

C 

EBBOLD-EBBOB 
EB20B=0. DO 
DO 46 1=1, H 

46 EBBOB=EBBOB*( (ZI (I)-Z(I))/ZE(I))**2 

F.PEOB=EBBOB/R 
PB IN 147, HR, EHBOl 

47 FOP.R AT (/* THE EBBOB FOHCTIOH BEFOBE ITEBATIOR MOflBEt*, 

1 12,* IS* ,712.8, 

IF (EBBOB. LT . EB BOLD) GO TO 49 

48 BEI=BEI-DP4 (1) 

BHH=BHH-D P4 (2) 
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tZ ET A= B2ETA -DP4 f3) 

OH K=K~DP4(4) 

1 GO TO 75 

2 49 CONTINUE ; 

C 

C SET OP AKD SOLVE THE SYSTEM OF LINEAB EQUATIONS. 

e 

3 DO 55 J=1,4 

4 B (J) =0. DO 

5 DO 50 1=1 , H 

6 50 B(J) =5 (J) ♦ (ZE(I)-Z(I)) *DZ(I,J)/ZE(I) **2 

7 DO 55 JJ=1,4 

8 A4 ( J , J J ) = 0. DO 

9 DO 55 1=1, B 


0 

1 



55 

A 4 (J,JJ) =A4 (J, JJ) ♦DZ (I , J) *DZ (I , JJ) /ZE(I) **2 
DO 58 J= 1, 4 

2 

3 

4 



58 

DP4 ( J) =1.000 
DO 58 JJ=1 , 4 

A 4 (J, J J) =A 4 (J,JJ)/B(J) /B (JJ) 

5 

S 

7 



59 

CALL DGELG (DP4,A4,4, 1, 1.E-14,IEB) 
DO 59 J= 1 , 4 

DP4 (J) =D P4 (JJ /B (J) 

3 

9 



60 

1 

PBINT6Q,IEB 

FOEMAT (/• THE ESBOB CODE FOB THE HATBII INVEBSION IS •, 

12,*.') 


o 

n n n 


ADJUST THE VALDES OF THE PABAHETEBS. 

) 

1 

l 




BEI=BEI*DP4 (1) 
B«N=BWN*DP4 (2) 
BZETA=BZETA*DP4 (3) 

i 


r 


K=K»DP4(4) 



c 


CHECK HHETHEB OB NOT ANOTHEB ITEBATION IS NECESSARY. 

\ 

j 


c 


DP(1) = DP4(1)/BZI 
DP (2) = DP4 (2)/B¥N 

1 

r 

1_ 




DP (3) =DP4 (3) /B ZETA 
DP(4) = DP4 (4)/K 
JJ=0 

l 

i 

i 



70 

DO 70 J=1, 4 

IF(DABS(DP(J)) .GT.1.D-3) JJ*1 
IF (JJ. EQ. 1) GO TO 71 

i 

i 



71 

PFIhT73,MM 
GO TO 75 

I F (H K. LT. 1 0 ) GO TO 34 

* 

• 



72 

1 

PBINT72, (DP (J) , J= 1 ,4) 

FOBMAT (/* 10 ITEBATIONS BAVE OCCOBED WITHOUT CON VEBGEBCE. •/ 

• THE PEBCENT CHANGES IN THE PARAMETERS ABE: ' 

2 //5I,4D14.5/) 

73 FOEMAT (/• CCBVEBGENCE OCCUBED OB ITEBATION hOMBEB',13,*. •) 

^ 74 FORMAT (*1«) 


c 

75 

CONTINUE 

PBIBT74 
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SDEI*BEI 

SBtiN-BHN 


SBZETA=BZEI 

SK=K 


PHASE 2 

DETERMINE INITIAL SET OP PARAMETERS 


BZETA=.05D0 
TZETA=SBZETA/2.D0 
TUN=BW N/2. DO 
IP(IBC) 129,236,130 

1 29 A=.25D0 

N E = 3 — — 

GO TO 132 

130 A=. 75D0 

NE=2 " — 

132 TMU=BEI/2.D0* {PI/BHN/BL**2) »*2* (SBZETA-BZETA) /A/TZET A»* (1.D0/NE) 
PRINT33 

PRI~NT7 7 ' 

PRINT78,BEI,BiN,BZETA,K 

PRINT79 

PRINT60, TrtU ,TWN,TZETA 

S AY E INIT I AL SET OF TI SSUE PA R AMETERS. 

STM0=TB0 

STUN=TWN 

STZETA =TZET A — 

THIS IS THE BEGINNING OF THE OUTSIDE LOOP. EACH RUN THROUGH THIS 

LOOP CONST I TUI IS ONE ilEEAIION. 

EBROR= 1. D20 
H M = 0 

134 KS=HB* 1 

CUECK EACH PABAflElllT FOR THE NON-N EGATT VIII CONDITION. 


_IP (TBU.LT. O.DO)GO TO 135 
I F ( TUN. L T . 0 .DO) GO TO 13 5 
IF (TZETA.LT. O.DO)GO TO 135 
GO TO 137 
MiiM = BB-1 
PRINT36,HBfl 
PRINT77 

PiiIhT7d, BEI ,BUN,BZ ETA, A 
PRINT79 

PRIKT80,TBO,TBN,TZETA 
TBU-STaU 
T BN- ST UN 
TZETA=STZETA 
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CALCPLATE Z AT EACH PBEQDERCT. 


BPH=BWS*2.D0*PI 

J PJM_= T *2 . D0»PI 

DO 145 1*1, II 

B3=DCHPLX (1.D0,2.D0*BZETA*P (I) /BPS) 

_TQ=_DCMrLI (I.DO^ 2.D0MZETA*P (I) /TPS) 

AEG= DC MPLX ( F ( I ) ♦Pl/TPS, Oj DO) /CDSQET (TQ) 

IP (IBC. EQ. 1 ) ARG=ABG/DCM PLX (2. DO, 0. DO) 

_COT=DCHPLX j 1. DO^O. DO) /CDTA K j[ A_RGJ 

LAttOA = CDSCBT(CDSQBl ( (DCB?LX*( ( (PI/BL) **2*P (I) /BPS) **2,0. DO) ♦ 
DCHPLX (P(I) ♦♦2*TflU/BEI,0. DO) /ABG/COT** IBC*I3C) /BQ) ) 

_ZA= L A H D A^DC ft PLXJ 8 A,0. DO) 

ZB=LAHDA*DCHPLX (EB,0.D0) 

ZL=LAHDA*OCRPLX(BL,O,D0) 

_Z 31 1 = DCft PLX (0- DO,-.5DO*P(I)/BEI) /LAB DA**3/BQ 
ZBI2=-CDSIKH (ZA) ♦COSISH (ZB) /CDSISH(ZL) 

—CDS I S (ZA) *CDSIN (ZB) /CDSIH (ZL) 

_ZBI= IUZB 12 

ZTI- DCNPLX (0. DO, P (I) /K) 

ZC=ZTI*ZBI 

_l -D R E A I. ( Z _C) 

Y=DlftAG(ZC) 

ZC= DCHPLX (1.D0,0. DO) /ZC 

Z (I) =C PARS ( ZC) 

CALCDLATE THE EEBIVITIVES OF Z AT EACH FBEQDEliCT. 


IP (1BC.EQ. 1 ) GO TO 
CSCS = '< CKPLX <1.00,0 

DZI(1)^DCKPLX <— P (I 

DZI ( 2) -DCHPLX (-12 J 
1 (DCHILX (1. DO, 0. DO) 

2 » (CSCS *C OT/ARG ) 

DZI (3)~sDCaPLX (0.D0 
1 * (COT + AEG*CSCS) 

GO TO 139 

138 CONTI SUE 
SECS-DCHPLX (1 . 00,0 

DZI (1) = DCHPLX_(l (I) 

DZI (2) "'DC ft PL I (-IKU 

1 (DCK PL X (1.DO,O.DG) 

2 • (SECS-DCHPLI ( 1.P0 
DZI (3) = DCK PLX (U.DJ 

1 * (Af:G*SECS-DCfi PLX ( 

139 COST1 S U£ 

i)Z I DL - ( DCM PLX (B A, 0 

2 ♦DCHPLX (EL, 0 

3 -DCHPLX (EL, 0 

5 A’CDCCSH ( 

DZIDL= DZIDL- (DCf.PL 
2 ♦DCHPL 


138 

.DO) /CD5IS (ABG) **2 

) ♦ *2/BEI/4.D0,0. DO)*COT/ABG/ 1 .AHDA**3 
«P (I) * ’ 2/TPN/B El/4 . DO", 0. D J ) /LA HD A * *3* 
-DCHPLX (0. DO, TPS ♦IZ ETA/P (I) /PI ♦♦2) *A8G**2) 

7-ihij*p (i) * iph 7pi^^27beT/4".ToT *Xe g/l a a d a • • T 


.DO) /CDCOS (ABG) **2 

♦♦2/BEI/4. DD,0.P0) /COT/A BG/L AH DA **3 
♦P (I) ♦♦2/1 PH/B Li/4 J DO, 0. D JJ/LAHDA**}* 

-DCS PLX ( 0. DO ,4.DQ*TPE # TZETA/P(I) /PI**2) *ARG**2) 
,O.DO) /CCT/ARG) 

, — 1 M U * P ( I) * T PE/P I ♦ *2/8 El ) ♦ARG/LAHDA**3 
1. DO, 0. DO) /COT) 

76" 0 f* C DCOf i'lifz a )* c'dsT SITfZ B ) r ' n 

.DO) •CDSIHIi(ZA)»CDCOSa (Z3) "7 L r. 

.DO) ♦CDSIKh(ZA) *>CDSISU (ZB) * 

ZL ) /C DJI hn (ZL) ) /CDS I Lu (ZL) 

X (BA , 0. DO) *00-03 (ZA) •CDSIH (ZB) 

X(E3,0. DO) *CDS1H (ZA) *CDCOS(ZB) 
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>3 


14 

5 

'6 


DP( 3) 5_PPJi 2XJS1 El A 

JJ-0 

DO 170 J=1 # 3 

170 IF(DA9S(DP(J) ) .GT.1.D-3) JJ=1 
IE (Jo. £Q. 1) GO TO 171 
PBINT73, MM 
GO TO 175 


171 IF (MM. LI. 10) GO 10 134 
PHINT172, (0 P (J ) , J= 1, 3) 

172 FOEMAK/ 1 1 0 I TERATIONS HAVE OCCORED BITHOUT CONVERGENCE. • i 

1 • THE PERCENT CHANGES IN THE PARAMETERS ABE: • 

2 //5X, 3D 14. 5/51,3*1 4. 5/) 

1 75 C O NTI N UE 

PBINT74 


PHASE 3 

DETEEHINE INITIAL SET OF PABAMETSBS. 

PRINT33 

PBINT77 

P R I NT 7 8 , BEI ,BBN , BZETA,K 

FBINT79 

PBINT60,TflU ,T8N ,TZE?A 


THIS IS THE BEGINNING OF THE CJTSIDE LOOP. 
LOOP CONSTITUTES ONE ITEBAT ION. 


EACH BUN THROUGH THIS 



ERROR= 1.D20 
MH=0 

234 BM=MM* 1 


CHECK EACH PARAMETER FOR THE HON *NEG AT I ¥ IT I CONDITION. 


IF (BEI.LT.O.DO)GO TO 235 
IP(BtfN.LT.O.bO) GO TO 235 
IF (TMU.LT. 0 . DO ) GO TO 235 
IF (T 4 N .LT. O.DO ) GO TO" 235 
IF (TZETA. LT.O. DO) GO TO 235 
IF (K. LI. 0. DC) GO TO 235 
GO TO 237 

235 MMM= MM-1 

PR INT36 , HMM 
PS IN 17 7 

FRINT78, BEI , BMH,BZETA,K 
PRIFT79 

PRINT8 0' # TMU,THN # IZETA 

236 BEI--SBEI 
BMN=SB*N 
BZETA- 3BZET A _ 

T A U = 0 . DO 
TWN^O.DO 

TZ ET A=~0. DO 
K=SK 

GO TO 275 


tJip W K) 
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237 CON TINUE 

C “ 

C CALCULATE Z AT EACH FEEQUEHCI. 

C 

BPN = BKN* 2 .D 0 *PI ~ " ~ 

TPN=T«N* 2 . DO*PI 

DO 2 45 1 = 1 . H • 

B 2 =DCMPLX ( 1 .D 0 , 2 .D 0 *BZETA*P(I) /BPN) 

TQ=DCHPLX( 1 .D 0 , 2 .D 0 * 1 ZETA*P (I) /TPN) 

ARG= DC MPLX (F(l) *PI/T PN , 0 - DO) /CDSQRT (TQ) 

IF (IBC.EQ. 1 ) ARG=ARG/DCMPLX ( 2 . DO, O.DO) 

CO T= DCHPLX ( 1 .D 0 , 0 . DO) /CDTAN (AEG) 

LAHOA=CDSQFT (CDSQRT ( (DCHPLX ( ((PI/BL) ** 2 »P (I) /BPN) *» 2 , 0 . DO) ♦ 

1 DCHPLX (P (I) **2 *IHU /BEI, 0 . DO) /ABG/COT*^ IBC*IBC) /BQ) ) 

ZA= LAS DA*DC HPLX (BA, O.DO) 

Z 3 = L AH DA *DC E P LX ( BB, O. DO) - 

ZL-L AB DA *DC HPL X (BL, 0 . DO) 

ZB 1 1 = DCM PLX ( 0 . D 0 ,-. 5 D 0 *P(I)/BEI) /LAHDA»* 3 /BQ 

Z 3 I 2 =CDSINH (Z A) ♦CDSINH (ZB) /CDSINH(ZL) 

1 -CDS IN (ZA) *CDSIN (ZB) /CDSIN (ZL) 

ZBI=ZBI 1 *ZB 12 

Z T I=DCM P LI ( 0 . DO ,P (I) /K) 

ZC-ZTI *ZBI 
I=DR EAL (ZC) 

Y=DIHA G (ZC) 

ZC=DCH PLX ( 1 . DO , 0 . DO) /ZC 
Z (I) =CDABS (ZC) 

C 

C CALCULATE THE DERIVITIVES OF Z AT EACH FREQUENCY. 

C 

IF (IBC- EQ. 1 )G 0 TO 238 

CSCS=DCMPLX (1 .DO, O.DO) /CDSIN (AEG) +*2 

DZI ( 1 ) =DCHPLI (iBU/q. DO* (P (I) /BEI)** 2 , 0 .D 0 ) *C 0 T/ABG/LAHDA **3 

DZI ( 2 ) =DCH PLX {-P<I)/EPN** 2 / 2 .D 0 ,O. DO) /LAHDA** 3 * 

1 (DCHPLX (P (I )/Bf h* (PI/BL) ** 4 , 0 '. DOj-DCHPLX (Q.DO,BZETA) ♦LAHDA** 4 ) 

DZI ( 3 ) =DCHPLX (-P (I) ** 2 /BEI/ 4 .D 0 , 0 . DO) *C 0 T/ABS/LAflDA **3 
DZI ( 4 ) =DCHPLX (-THO*? (I ) ** 2 /TPN/ 3 EI /4 . DO, 0 . DO ) /LA HDA** 3 * 

1 (DCM PLX ( 1 . DO, 0 . EQ) -DCM PLX ( 0 . DO , TPN *iZET A/P (I) /PI** 27 *ARG** 2 ) 

2 * (CSCS+COT/ARG) 

DZI ( 5 ) =DCMPLX (O.DO,-TBU*P(I) *TPN/PI** 2 /BEI/ 4 . DO) * ARG/L ABDA **3 
1 *(C 0 T*A 2 G*C SCS) 

GO TO 239 

238 C OKT INUE 

SEC 5 =DCflPLX (1 .DO, O.DO) /CDCOS (ARG) **2 

DZI ( 1 ) =DCf!PLX (-XHD/ 4 .D 0 * (P (I ) /BEL) ** 2 , 0 . DO) /C 0 T/ARG/LASDA **3 
DZI ( 2 ) =DCBPLX (-P (I)/£PN** 2 / 2 . DO, 0 . DO) /LAMDA** 3 * 

1 (DCHPLX (? (I J/BPa* (PI/BL) ** 470 . DO) -DCS PLX (OTdO, BZEIA) •LAB'DT*^! 

DZI ( 3 ) =DCMPLX ( P (I) ♦♦ 2 /BEI/ 4 . DO , 0 .D 0 ) /C 0 I/ADG/LAHDA **3 
DZI ( 4 ) =DCKPLX (-TBU *P (I ) ** 2 /TPN/F> EI/ 4 . DO , 0 . DO ) /LABDA** 3 * 

1 (DC H P Lll 17 D OTOTC 0 ) ^ITC M r L X" ( 0 .’DO , 2 } ^ D 0 *-^P N"AT 2 E I K 72 JIT 7 P 1 * * 2 T*ARG 

2 * (SECS-DCBPLX ( 1 . DO, O.DO) /COT/AEG) 

DZI ( 5 ) -DCHPLX (O.DO,-IMU*P (I) *IPN/PI** 2 /BEI) *I.BG/LAHDA **3 

1 * ( A 2 G-? SECS- DCHPLX ( 1 . DO D 0 ) /COT) — — 

—239 CONTINUE 

DZID 1 = (DCBPLX (BA, O.DO) »CDC 0 SB(Z 1 ) »CDSIMH (ZB) 




i 




AMIft irlfhriai irfc 


co cr> <n nl -r in -n!r^ co 
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♦DCflPlXJBB^O.PO) *CDSJNH(JA).*f:nCOSHJZB) 

-DCftPLI (BL,O.DO) ♦CDSIN:-;(ZA)*CDSINH (ZB) 

* CDCOSH (ZL)/CDS1NH (ZL) )/CDSlNU (ZL) 

DZIDL= DZIDL- ( DCMPLX (BA,O.D0) *CDCOS (ZA) *CDS1N (ZB) 
♦DCKPLX(BB,0. DO) *CDSIN (ZA) *CDCOS (ZB) 

- DCMPLX (BL, 0. DO) *CDSIN (ZA) *CDSIN (ZB) 

*CPCOS (ZL)yCDSIN (ZL) ) /CDSIN (ZL) 

DZIDL=DZIDL*ZBI HDCMFLX (-3. DO, 0. DO) /LASDA*ZBI 

DO 240 J=1 , 5 

DZI(J ) - DZI (J) /BQ*D ZI DL 

DZI (1) =DZI ( 1) -ZEI/DCKPLX (BEI,0. DO) 

DZI (2) =DZI (2) ♦ ZEI/BQ* DCMPLX (0. DO, 2. D0*BZETA*P (I) /BPN**2) 
DZI (6) =-ZTI/DCMFLX (K,O.DO) 

DZI (2) =DZI (2) ♦DCMPLX (2.D0*PI ,0. DO) 

DZI ( 4) =DZI (4) ♦DCMPLX (2. D0^PI,0. DO) 

DO 245 J=1,6 

DX (J) =DK EAL (DZI (J) ) " 

DY (J) = DIMAG (DZI (J) ) 

DZ (I , J) =-2 ( I) *»3» (DX (J) *X+DY (J) *T) 

CALCULATE AND PBI NT THE EBEOB FUNCTION. 

~ERBOLD=EREOi 
EEBOB=0. DO 
DO 246 I = 1,N 

£EROR= EBKQB ♦ ( (Z£ (I) - Z (I) ) /ZE (I) ) **2 
ERROE= EBHOE/N 
PR I NT 4 7, MM, ERROE 
" I F( ERE OR. LT^ESEOLD ) GO TO 249 
B2I=BE1-DP6 (1) 

DW N- BW N-DP6 (2) 

TM U= TMU-DP6 (3) 

TifN = TiiN-DP6 (4) 

TZETA~TZETA-DP6 (5) 

X=K— DPb (6) ’ 

GO TO 275 

CONTINUE 


SET DP AND SOLVE THE SYSTEM OF LINEAB EQUATIONS. 


-Q 


DO 255 J*1,6 
B (J) =0. DO 
DO 250 1=1, N 

B ( J) =b (J) ♦ (ZE (I)-Z (I)) ♦DZ{I # jjVZE(I) **2 
DO 255 J J= 1 ,6 
A6 (J,JJ) =0. DO 
DO 255 I=T, N 

£6 (J, JJ) =A6 (J, JJ) ♦DZ (I , J) *DZ (I , J J) /ZE(I) *»2 
DO 258 J = 1 , 6 
DP6 ( J) =1TDD0 
DO 258 JJ= 1 ,6 

A6 (J^JJ) =A 6 (J,JJ)/B(J)/B(JJ) 

CALL DC ELG (D?6,A6,6, ?717E- 1 47IEB] 

DO 259 J=1 , 6 
DPS (J) -DPS (J)/B(J) 
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PRINT6 0. IE B 

C ~ 

C ADJUST THE VALUES OF THE PARAMETERS. 

C 

BEI=BEI*DP6 (1) ~ " 

BV N= BW N* DP6 (2) 

TH tl= I M tUD P6 (3) 

T«N=TWN'DP6 (4) 

TZETA = TZETA *DP6 (5) 

K=K+DP6 (6) 

C 

C CHECK HHETHEB OB NOT ANOTHEB ITERATION IS NECESSABT. 

C 

DP(1) = DP6{1)/BEI 
DP (2) = DP6 (2) /BHN 

DP ( 3) = DP 6 ( 3 ) / T MO 

DP(4) = DP6(4)/TBN 
DP (5) = DP6 (5)/TZETA 

DPJ 6) -DP6 (6)/K 

JJ = 0 

DO 270 J=1 , 6 

270 IF (D A B SJ DPJ J) ) . GT. 1. D-3) JJ-1 

IF (J J. EQ. 1) GO TO 271 

PRINT7 3, Bfl 

GO_TO J75 

27 1 IF (MM. LT. 10) GO TO 234 
PRINT 1 72, {DP (J) ,J=1,6) 

-275 CONTINUE 
PRINT74 
C 

C PRINT THE FINAL PARAMETER VALUES. 

C 

PRINT76, TITLE, EL, DRATIO 

76 FORMAT (5X/1 5A4//5X , • BON E LENGTH ' ,5 1, * PROBE LOCAT ION '/F9. 1 , Fl 5. 1/) 
PRINT? 7' 

PRINT? 8, BEI , BVN,BZETA,K 
PR I NT 7 9 

PS IN 'Id 0 ,TMU,T BN , T Z E T A 

77 FORMAT (51, • EON E STIFFNESS* ,5X,' BOSE EAT FBEQ' # 5X, 

1 'BONE DAMPING', 5X, 'SKIN STIPFNESS') 

78 FORM AT (D 16. 5^ 1 3*1 ,F19. 4,022. 5/) 

79 FORMAT (5X, 'TISSUE DASS/LENGTH* ,51, 'TISSUE NAT FBEQ',5X, 

1 'TISSUE DAMPING') 

80 FOB MAT (F 10. 2, F23. 1 , F21. */) 

PRINT 81 

81 FORMAT (3 IX, 'EXPERIMENTAL', 14X,« THEORETICAL' /4 IX, 'PHASE' ,211, 

1 '"P HA'S E V 13x7* ? ELQ T 7B‘X 7 1 1 a? LD A'NC A N G LE t , 7 X , • I iTPED AaCE^ 

2 51, 'ANGLE'/) 

C 

C RECAlC UIXT E~T H r~lM 'PEDTSITE i — 

* ■** *'*'.*, **_ 

C ir 7 • * iK fT", 

B? N= BV N* 2. D 0*PI 

tpn=tnn*T.dO*p! 

DO 85 1=1, i 

BQ=DCM PLX ( 1 . DO , 2. D0*BZETA»P (I) /BPM) 
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C 

c 
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JPiTH’J.EQ.O.DO) GO TO_82 

T0=DCHPLX ( 1 .DO ,2. DO* TZETA* P (I) /TPN) 

AR G= DCMPL X ( P (I ) *P1/TPK,0.D0)/CDSQRT(TQ) 

IF(IBC. EQ. 1) APG^= ARG/CCM PLX (2. D0,0. DO) 

COT= DC MPLX (1. D0,0.~ CO) /C DT AN (AEG) 

LAHDA = CDSQ8T(CDSQB1 ( (DCBPLX ( ((PI/BL) **2*P(I) /BPN) **2,0. DO) ♦ 

1 D C M P L X (P ( I) **2*IBU/ B EI , 0. D O) /AB G/C OT*»IBC*I BC) /BQ) ) 

GO TO 83 

82 LABDA=DCBPLX (PI/BL, 0. DO) *CDSQET (DCMPLX (P (I) /BPN, O.DO) /CDSQFT (BQ) ) 

83 ZA =LA MDA*DCBPLI(BA,0 .DQ) 

ZB=LAM DA *DC flPLX (BB ,O.DO) 

ZL=LAMDA*DCBPLX (BL,O.DO) 

Z9I=DCBPLX {0. DO,-. 5D0*P (I) /DEI) /LAMD A**3/B0* 

1 (CDSINH(ZA) *CDSINH (ZB) /CDS'lNU (ZL) -CDSIN (ZA) *CDSIN (ZB) /CDSIN (ZL) ) 
ZTI=DCBPLX (O.DO,P(I)/K) 

ZO DCM FLX (1 - D0,0. DO) /(ZTI*ZBI) 

Z ( I j =C DA DS ( ZC ) 

PBI (I) =DATA N (DIBAG (ZC) /DBEAL(ZC) ) *180. DO/PI 


PRINT THE IBPEDANCE. 

85 PRINT86, I, N (I)jZP(I) ,PHIE^I) ,Z (I) ,PHI (I) 

8b FORMAT (17, F 13. 2, D1 6. 4, F10. 2, D16. 4, Ft 0. 2j 

CALCULATE AND PR IK T THE EE JR O R FUNCTION. 

ERROR - 0. DO 

P0_87_I=1,N 

87 _ EBRC»R = ERROR* ( (ZE (I) -Z (I) ) /ZE (I) ) **2 
EHhOB= EEROR/N 
PBINT88, ERROR 

83 FOE MAT (/• THE ERROR FUNCTION FOB THIS SET OF PARAHETEBS IS^ 

1 F12.8, •. »/) 

PRINT 74 


PLOT THE IB PEDANCE • 


DO 90 1=1, H 
NP(I) = SN GL ( W (I ) ) 

ZEP (I) =SEGL (ZE (I) ) 

PHIEP { I) — S N GL (f HIE (I) ) 

ZP (I) = SNGL (Z (I) ) 

90 PHIP (I)=SNGL(PHI (I)) 

CALL P LTOF S (1 . , 1./2. ,3.,1./2.,1.5,4.5) 

CALL PLGAXS (1. 5,4. 5, • FBEQU EKCT* ,-9,6. ,0. , 1. , 1. /2. ) - 
CALL PLGAXS (1. 5,4. 5, * IM PED AN CE » , 9, 6. , 90. , 3. , 1./2. ) 

CALL PGR ID ( 1. 5,4.5, 2., 2. ,3,3) 

CALL PLTLOG (3) 

CALL PLIN£(1TP(1) ,Z P( 1) , H, 1 , 0 ,0, 1 ) 

CALL PL1EE («,: (t) ,ZE? (1) ,N, 1,-1,“07TJ 
CALL PLTBEC 

CALL PLTOFS (1. ,1./2. ,-9 0. ,90. , 1. 5, 1. 5) 
call plgaxs'(V. 5, l.37» iEEQUEECi 1 ,-97^'T7^7TrrT.72t7 
CALL PAX1S (1. 5,1. 5, 'PHASE AEGIE’,11,2.,90. ,-90.,90.,„25) 
CALL PGEID (1.5,1.5,2.,1.,3,2) 
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! CALL_PI.TLOG (2_) 

t V ,* CALL PLINZ(feP (1) ,PKIP(1) ,H, 1,0, 0,1) 

« ' CALL PLINE(VP(1),PHIE?(1),S, 1,-1, 0,1) 

■ . CALL PLTBEC 

CALL PGEID(0.,0.,8.5,11.,1,1) " 

r - CALL PSYHB(1. 5,. 5,. 125, TITLE (1) , 0. ,60,0) 

C A LL PL TEHD 

* END T ~~ 

■ TIONS IN EFFECT* ID, EBCDI C, SOUBCE, NOLIST, NODECK, LOAD, NOflAP 

1 riOKS IS EFFECT* HARE = RAIN , LINECHT = 57 

\ All STICS* SOOECE STATEHENTS = 5Q9,P£0GHAf! SIZE = 35900 

f ATI STICS* BO DIAGNOSTICS GENEBATED 

i ES IN RAIN 
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EAL.XONC.TI ON _PPEAL*0j 
COMPLEX* 16 i',DCfiPLX 
RL AL*3 Y,CDADS,EBLI 
DR EAL=CDAB5 (( X + DCOKJG 
I = DBLE (R EAL(X) ) 
D8EAL=DSIGH (DREAL, T) 
RETURN 


[2.D0-0. DO] 


PIICHS IN EFFECT* I D, EBCDI C, SO OECE, NOLIST, NODECK, LOAD, BOM A P 
PTIONS IN EFFECT* KANE - DREAL , LI NEC NT *= 57 

TATI STICS* SOURCE STATEMENTS = 8,PB0'gRAH SIZE = 

I A7ISTICS* HO DIAGNOSTICS GENERATED 
OES IN DBEAL 




r 127 -p 

\ GAN TESBINAL SISTEM FOBTEAN G(91336) DIB 1G 09-19-78 2 


1 1 

BEAL FUNCTION DIMAG* 8(1) 


2 l ,.y 

3 
• U 

COMPLEX* U» X , DCKPL X 
8EAL*8 I , CD ADS , DELE 

DIHAG=CDABS ( (X-DCOXJG(X) ) /DCSPLX (2.DO,O.DO)) 


5 

I-DBLE (AIKAG(X) ) 


6 

DIBAG-DSIGN (DIBAG, I) 


7 

BETORN 


8 

END 



PIIONS IN EFFECT* ID, EBCDIC, SOURCE, NOLI ST, NODECK, LOAD, NOB A P 

?T 1 0 N 5 _ IN EFFEC T* NAM E - DIH AG , LIHECHT = 57 

TATI STI C S* SOUBCE' STATEMENTS = 6 , PBOGEAB SIZE * 

TATISIICS* SO DIAGNOSTICS GEHEEATED 
Or.S IN DIMAG 


127-Q 


AN TEESIN AL SISJEB FOBTBAN G (41336) 


CDSIVH 


09-19-78 


L r > CO 3PlEX.._r U NCTXUN CDS INH »JUX) 

COMPLEX* 16 X,CDEXP,DCiJPLX 

i CDSINU= (CDEXP (X) -CEEXP (-X) )/DCHPLX (2-DO,O.DO) 

BETUBN 


* END 

’HONS IN EFFECT* ID, EBCDIC, SO DEC E, NOLI ST, NODECK, LOAD, B0H1P 

>n ON S IN EFFECT* N AM E = CCS INH . LIMECST ■ 57 

'AIISIICS* SOUECE STATEBEHTS » .5, PBOGBAH SIZE « 

iTISTICS* NO DIAGNOSTICS GENEBATED ^ 

)ES IN CDSINB 
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IAN TEBHINAL STSTEfl FOBTBAH G (41336) CDCOSB 09-19-78 21 

CO.fi P.l l X_IJ) H CTIOH CDCOSH »16m 

:U COMPLEX* 16 X,CDEXP,DCHPLX 

I CDCOSH= (CDEXP(X) *CDEXP (-X) ) /DCHPLX (2 • DO* 0* DO) 

. .RETOBN 

► END ‘ 

>TI ONS IN EFFECT* ID, EBCDIC, SOURCE, NOLIST, NODECK, LOAD, MOB AP 
>TIONS IN EFFECT* NABE * CDCOSH . LINECNT * 57 

'AIISTICS* SODECE STATEMENTS = 5, PROGRAM SIZE * 522 

ATISTICS* 80 DIAGNOSTICS GENEBATED 
)ES IN CDCOSH 


GIN TERMINAL STSTEH FOBTBAH 6(41336) 


COTAS 


127-S H 
09-19-78 2. 


1 r*, XOBPJL£X._PJ] N CTXON _CJDT A N *_16_(I) 

2'«- / COMPLEX* 16 X,CDSIN,CDCOS,DCHPLX 

3 IF (DIM AG (X) .LE. 34. Cl ) GO TO 1 

4 ' CDTAN* DCMPLX (0. DO, 1. DO) 

5 CO TO 3 

6 1 IP(DIMAG(I) .CT.-34.D1)GO TO 2 

7 CDTA N= DC B PL X (0. DO. -1 » DO) 

3 GO TO 3 

9 2 CDTAN*CDSIN (X) /COCOS (X) 

0 3 BBTOBM 

1 END 

PTION5 IN EFFECT* ID, EBCDIC, SOURCE, NOLIST, NODECK, LOAD, MOUAP 

PI 10 NS IN EFFECT* G A W E = CDTAS , LI NEC NT « 57 

TATISIICS* SDUBCE STATEUENTS ■ 11,PB0G£AH SIZE * 

TAIISTICS* VO DIAGNOSTICS GENEBATID 
3RS IN CDTAM 
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&. I£S ULTS OF I£ 7170 TESTS 08 THE F0RE1RHS OP SEFEW HO HIM 
SUBJECTS 

The results of the in vivo tests perforned on Subject TT 
are presented and discussed in Section YII.l. Siailar results 
froa seven other subjects hare been obtained and are presented 
here. Driving-point aecbanical iapedance plots associated with 
•10C< 500 and 600 gran-force preloads are given for each subject* 
The corresponding paraaetric values in each case are listed in 
Table 7.1. 


o 
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TlBtES 



TABLE 2.1 


Three Basic Types of Mechanical Eleients 



■ass 

daup >r 

spring 

Equation of Motion 

f * mi 

f = cx 

f « kx 

F/a 

a 

C/P 

*/P* 

Slope on log-log plot 

00 

-450 

-63.40 

P/t (iopedance) 

»P 

c 

*/P 

Slope on log-log plot 

1*50 

0® 

-450 

F/6 

■P 2 

cp 

k 

Slope on log-log plot 

63.4® 

450 

OO 


TABLE 3.1 
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Boundary Conditions 


at x * 0 


1. Simply-supported 


Ti * 0 
H, * 0 


It * 0 


2. Botational spring 
on one end 


T. * 0 

k, 0, - H, 


lz * 0 


a 2 * 0 


3. Botational spring 
on each end 


y, * 0 


k, 0, - B, 


0. Translational spring 
on one end 


k, 1 , ♦ *. 

a, * o 


lz * 0 


a 2 = o 


5. Translational spring 
on each end 


k, J, ♦ ▼, 

a, *= o 


a, * o 


at i * -e 


6. Translational spring 
on an extended beam 


k ils ♦ *3 


y* * 0 


O 
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TABLE 3.2 

Non-diacnsional Parameter Definitions 


Mon-diaensional 

Paraaeter 


Definition in teras 
of Model Paraaeters 


5 y * ujYl/22 

a a/L 

p b/L 

H p^/p = L*/jt* p^U)*/EI 

B o/lt'f 

S k/K = kL 5 /48EI »* 

T 2kL»/EI ** 

B AL/2EI >* 

£ 

C y CfUj/k 

C* c^uj/k^ 


» *k is the spring constant of the spring in series with the 
bean. 


**k is the spring constant of the translational spring at a 
support. 


* * k is the spring constant of the rotational spring at a 

support. 
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TABLE 4.2 

Static Stiffnesses for Beaas iith Various Boundary Conditions 


listed belov for several different 


The stiffness of a bean is 


K * ♦ 3EIL/a*b* 


vbere expressions for 4 are 
boundary conditions. 


Boundary conditions 


1. Sinply-supported 


2. Rotational spring 
on one end 


3. Rotational spring 
on each end 


4. Translational spring 
on one end 


5. Translational spring 
on each end 


6. Translational spring 
on an extended bean 


♦ 

1 

6+4B, 

6*«{3a*40) R, 

6*4R, ♦4R 2 *2R, R t 

6+tt(3«*4£)B, * 0(3£*4a) B l *2a^B l B z 

T,«* 

6*T,x* 

T, T 2 **p* 

6(T,a* + T 2 p2)*T t T 2 «*|J* 

24+4£»T 3 *4e*T s 

24t4£>T 3 *A(3a*4p)£2T 3 


135 


TABLE 6.1 

Paraaetric Values for the Foreara of Honkey 663 


Paraaeter 

laae 

Olnar report length 
Length-to-probe-location ratio 
Olnar bending stiffness 
Ulnar fnndaaental frequency 
Olnar daaping ratio 
Support rotational stiffness 
Support rotational daaping 
Tissue aass per unit length 
Tissue fnndaaental frequency 
Tissue daaping ratio 
Skin stiffness 

Condition 
Excised ulna 
Busculature reioved 
Probe on ulna 


Paraaetric 

Syabol Value 

I 17.1 ca 

a 0. 6 

El 2.9795x10* dyne ca 2 

Co 332.0 Bz 

y 0.0425 

k, 0.86535x10* dyne ca 
c, 1.7136x10 s dyne ca s 

Pc 1-85 g/ca 

w* 174.0 Hz 

0.4050 

k 2.2098x10* dyne/ca 


Value of Error Punction 
0.0086 
0.0112 
0.0115 


Intact ara 


0.0132 


136 


TABLE 6.2 

Paraaetric Values foe the Foreara of Bonkey 665 


Paraaeter 

Ease 

Ulnar support length 
Length-to-probe-location ratio 
Ulnar bending stiffness 
Ulnar fondaaental freguencj 
Ulnar daaping ratio 
Support rotational stiffness 
Support rotational damping 
Tissue aass per unit length 
Tissue fundamental frequency 
Tissue daaping ratio 
Skin stiffness 

Condition 
Excised ulna 
Husculature removed 
Probe on ulna 


Paraaetric 

Syabol Value 

L 17.2 cb 

ft 0.6 

El 5.0311x10* dyne ca* 

to 350. Q fix 

* 0.0364 

k, 4.4382x10* dyne ca 

c, 2.2225x10 s dyne ca s 

6.96 g/ca 
tOj. 101.0 B z 

0.0792 

k 1.1155x10* dyne/ca 


Value of Error Function 
0.0226 
0.0127 
0.0208 


Intact at a 


0.0653 
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T1BLE 6.3 


Pare™etric Values for the 
Parameter 

Mate 

Dinar support length 
Length-to-probe-location ratio 
Dinar bending stiffness (ED) 

Dinar bending stiffness (HH) 

Dinar fundaaental frequency 
Dinar daaping ratio 
Support rotational stiffness 
Support rotational daaping 
Tissue aass per unit length 
Tissue fundaaental frequency 
Tissue daaping ratio 
Shin stif fness-400 gn preload 
Skin stiffness-600 ga preload 


Foreara of Donkey 659 
* 

Paraaetric 

Syabol Value 


L 

17.2 ca 


Oc 

0.6 


ZI 

5.2498x10* 

dyne ca? 

El 

7.7120x10* 

dyne ca 2 

<0 

377.4 Hz 


r 

0.0267 



3.4682x10* 

dyne ca 

C, 

4.2095x10* 

dyne ca s 

p* 

4.02 g/ca 



145.1 Hz 



0.5827 


k 

1.3543x10* 

dyne/c a 

k 

1.3924x10* 

dyne/ca 


Condition 

Excised ulna 

Husculature reioved 

Probe on ulna 

Intact ara 400 ga preload 

Intact ara 600 ga preload 


Value of Error Function 
0.0094 
0.0179 
0.0111 
0.0122 
0.0169 


TABLE 6.0 


Bending Stiffness Heasureaents on the Ulna of Bonkey 659 


Test 


Bending Stiffness 
El (10* dyne ca*) 


dpbi test (auscolature reaoved) 7.712 
DPBI test (excised ulna) 5.246 
Percent difference 32. OX 
Three-point bending test (HTS aachine} 4.827 
Percent difference 37. 4% 
Bepeat bending test on dry bone 4.530 


Percent difference 


41. 3X 


TABLE 6.5 


Hechanical Properties of the Aluninon Bean 


Bending Fundaaental 
Stiffness Frequency 

El 40 


(10* dyne cn*) (Bz) 

Predicted values 5.587 429.2 

Corrected for enlarged ends 5.670 448.2 

Percent difference 1.5% 4.4% 

Measured values 6.090 489.3 


percent difference 


9.0% 


14.0% 


Beasared Quantities Bone Paraoeters Tissue Paraaeters Skin Parameters Calculated 


Bone Bass/length 

CM 

in 

CM 

H (g/ca) 

• 

r- 

• 


r* 

CM 

Stiffness 600 ga pi 
k (10 s dyne/ca) 

Cl 

ci 

a 

0 

0 

vO 
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m 

in 


CM 
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CO 

Stiffness 500 ga pi 


CO 
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X (10 s dyne/ca) 
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Stiffness 400 ga pi 

ip 

co 
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IP 

k (10 s dyne/ca) 
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CM 
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CT> 

Damping ratio 

w> 

0 

(fj 

o 

r* 

CO 

• 

o 

(diaensionless) 

• 

O 

Fundamental frequency 
u > f (Hz) 

137.8 

162.9 

Hass per unit length 
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Pf (ga/ca) 
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CO 
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Fundamental frequency 
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Figure 1.1. Buaan Long Bones. 


(a) ira and foreara shoving relative size, shape and position of 
its bones, (b) Thigh and leg shoving relative size, shape and 
position of its bones. 
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Figure 1.4. Saaple Output Proa Thompson's Prograa. 
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Figure 2.1. Orae’s First flodel of the Ulta in Thompson^ 
Experiaental Procedure. 
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Figure 3.1. Diagrams of Beaa Hodels 


(a) Case 1: siaply-snpported 

(b) Case 2: rotational spring on one end 

(c) Case 3: rotational spring on each end 

(d) Case 4: translational spring on one end 

(e) Case 5: translational spring on each end 

(f) Case 6: translational spring on an extended bean 





Pi gore 3.3. The Foundation. 


The coordinate systea and boundary conditions: (a) fixed, (b) 
free. 
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figure 4.5. DPMI of Case 2: Botational Spring on One End 
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Figure 4.7. DPHI of Case 4: Translational Spring on One End 
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Figure 4.10. DPMI of Cases 
diaeosionalized. 
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Figure 4.12. Taper. 
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rigore 4.13. DPHI Plot Exhibiting the Dependaoce of the 
lass Pec Onit Length of a Fixed Foundation. 


(The bean boundary conditions are siaply-aopported 
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Figure 4.14. DPHI Plot Exhibiting the Dependence of the 
Daaping Batio of a Fixed Foundation. 


(The beaa boundary conditions are siapl y-supported.) 







\ Figure 4.15. DPHI Plot Exhibiting the Dependance of the 

t Hass Per Unit Length of a Free Foundation. 

(, (The tea* boundary conditions are siiply-supported.) 
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C Figure 4. 17. Conparison Between Actual Biniaun DPMI and 

Approxiaate Equations. 


(a) Equation (4.17), fixed foundation, (b) equation (4.18), free 
foundation. 
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Figure 5.2. Flow Chart of the Conputer Pro gran 
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Figure 5.5. Saaple Output Froa coaputer Prograa. 
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Figure 6.1. Bonkey lea ia Test Fixture. 


(a) latact era, (b) Probe on ulna, (c) Husculatore reaoved 
Excised ulna. 





















































Figure 6.16. Bending Fixture Used For Three-point Bending 
Test on the Ulna of Honkey 659. 
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Figure 6. 19. Dimensions of the Aluminas Beam and its 
Support Brackets. 
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Pigare 6.20. DPBI of the Aluainua Beai 
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Figore 7.5. DPHI of HonJcey 16: Tibia. 
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Figure B.l. The Fleaents of a Tapered Bean 


(a) Linear taper, (b) Quadratic taper 
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Figure C.l. True Biniiua of a Discrete DPBI Plot. 
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